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SUMMARY 


The NASCAP ( NASA Charging Analyzer Program) code 
simulates the charging process for a complex object in 
either tenuous plasma (geosynchronous orbit) or ground test 
(electron gun source) environment. NASCAP is written in 
FORTRAN V and runs under EXEC-8 for UNIVAC 1100 series com- 
puters. This volume contains detailed specifications 
needed to run or modify the code. The physical content of 
NASCAP is described in NASA CR-135258. 

The object definition section, OBJDEF, allows the 
test object to be easily defined in the cubic mesh. The 
test object is composed of conducting sections which may 
be wholly or partially covered with thin dielectric coatings. 
It is built of rectangular parallelepipeds and objects 
having slanted surfaces (WEDGE, TETRAHEDRON) . Complex 
objects such as an octagonal right cylinder (OCTAGON) and 
a twenty-six faceted quasi-sphere (QSPHERE) may be ex- 
plicitly defined. The different conducting sections may 
be floating, held at fixed potentials, or biased relative 
to one another. 

The potential section, POTENT, obtains the electro- 
static potential in the space surrounding the object. It 
uses the conjugate gradient method to solve the finire 
element formulation of Poisson's equation. The potential 
is defined in the inner mesh containing the test object, 
and in an arbitrary number of nested outer meshes having 
successively doubled grid spacings. (This enables treat- 
ment of a large volume of space exterior to the object.) 

To promote efficiency of the potential solver at each time- 
step, several potential scaling options are available to 
produce "initial guess" potentials based on a previous set. 



The CHARGE section of NASCAP treats charge re- 
distribution among the surface cells of the object as 
well as charging through radiation bombardment. In the 
ground test case the incident flux is determined by 
forward particle tracking from an electron gun of speci- 
fied beam characteristics. In the plasma case the inci- 
dent flux may be determined by the reverse trajectory 
sampling method, or by a Maxwell probe approximation. 

In all cases energy- and angle-dependent formulations for 
electron backscatter and secondary electron emission due 
- to incident protons and electrons are used. A full 
shadowing treatment is used in computing photoemission 
due to incident sunlight. A first-order explicit photo- 
sheath treatment is available, and a "LONGTIME STEP" 
option allows timesteps of one minute or more during dif- 
ferential charging. 

NASCAP has facilities for extensive graphical out- 
put, including several types of object display plots, 
potential contour plots, space charge density contour 
plots, current density plots, and particle trajectory 
plots. 
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1 . INTRODUCTION 


1.1 WHAT IS NASCAP? 

The NASCAP code was developed by Systems, Science and 
Software under contract to NASA-Lewis Research Center to sinu 
late the charging process for a complex object in either tenu 
ous plasma (geosynchronous orbit) or test tank (electron gun 
source) environment. A block diagram of the NASCAP code is 
shown in Figure 1.1. The focus of the code is the CHARGE 
section, which calculates the fluxes to the surface elements 
of the object and the resulting material response. The 
particle tracking capabilities of the code are included in 
the CHARGE section. The POTENT segment is a conjugate gradi- 
ent (CG) Poisson solver for the electrostatic potential on 
and about the object. The potential may be calculated in an 
arbitrarily large volume through the use of successively 
coarser outer grids. The OBJDEF segment allows definition of 
an object composed of bare or dielectric-coated conductors. 
The surfaces of the object may be normal to any of the twenty 
six symmetry directions of the cubic mesh. The IIIDCEL seg- 
ment calculates the portion of each surface cell exposed to 
sunlight. Graphical capabilities of the NASCAP code include 
object display, particle trajectories, and contours of 
electrostatic potential, space charge, or test tank current. 

NASCAP has been written in FORTRAN V and runs urder 

EXEC-8 system for UNIVAC 1100 series computers. It has been 

3 ^ 

tested and installed on the UNIVAC 1108 at S using the S' 
graphics library, and on the UNIVAC 1110 at NASA-LeRC using 
both the FR80LIB and LEWIGS graphics libraries. Core storage 
requirements are minimized by extensive segmentation and use 
of fast block transfer between core and scratch files. 

NASCAP runs in approximately 7 OK words (decimal) of core. 

We believe NASCAP to be transferable to any medium to large 
size computer having a word length of 36 bits or longer. 
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Figure 1.1. Block diagram of the NASCAP code 
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1.2 WHAT DOES NASCAP DO? 

NASCAP simulates the charging process in an explicit 
time-stepping fashion. Suppose that at a particular tine we 
know the charge distribution on the object, the electrostatic 
and magnetostatic fields about it, and a description of the 
environment. NASCAP first calculates the flux of electrons 
and ions incident upon each surface cell of the object. In 
the test tank case this is done by tracking particles forward 
from the electron gun source, whose characteristics are known 
from calibration data. In space this may be done either by 
approximating each surface element as a surface element of 
a spherical probe in an isotropic Maxwellian plasma, or by 
sampling the plasma distribution using reverse trajectories 
frjm the satellite. The latter method automatically takes 
particle shadowing into account. Next, electron backscatter 
and low energy electron emission are taken into account, 
including a full treatment of sunlight shadowing. Optionally, 
currents and space charge due to the low energy emitted 
electrons may be calculated by particle tracking. 

The fluxes calculated above are assumed to hold for a 
user specified time, and the charges on the satellite are up- 
dated accordingly. The second portion of a time step is the 
calculation of the electrostatic field about the object in 
the presence of the new charge distribution. The solution of 
Poisson's equation about a complex three-dimensional object 
is a large task which can easily come to dominate the cost of 
a charging simulation. NASCAP uses a sophisticated conjugate 
gradient iterative scheme to solve a finite element formula- 
tion of Poisson's equation. Startup, restart, and potential 
scaling features enable the user to provide a good starting 
point for the iterative process, and thus minimise the number 
of iterations required to achieve a potential field with 
tolerable levels of error. When the new potential has been 
calculated, the charging process for the next time step can 
begin . 
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1.3 


WHAT DOES NACCAP EXPECT FROM ME? 


Because NASCAP is a very general, very flexible code, 
a great deal of user input is required: 

1. The object to be tested must be defined geo- 
metrically and electrically. 

2. The properties of the surface materials of the 
object must be characterised. 

3. The environment in which the simulation is to 
take place must be defined. 

4. The necessary computational resources must be 
made available. 

5. The user must determine the timescale of the 
simulation. 

6. The user must provide parameters concerning 

the size of the computational space and the level 
of accuracy of various parts of the simulation. 

7. Intelligent choices must be made among the 
various options provided. 

3.4 HOW DO I USE NASCAP? 

Because NASCAP consists of several complex modules 
combined into one very complex code, it is not easy to use. 

In this section we discuss some strategies and tactics for 
carrying out a NASCAP simulation, postponing detailed input 
specifications to subsequent chapters. 

1.4.1 Computational Space and Object Definition 

NASCAP focuses on an "inner grid" of dimension 
17 x 17 x 4n+l (4 £ n <_ 8) within which the object is .located. 
Outer grids have similar dimension but successively double 
the physical mesh spacing. The first task of the user is to 
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determine the smallest computational space which allows for 
adequate representation of the object, and thus minimize -the 
effort and expense required for the simulation. A simple 
object can often be run as a one-grid problem. A typical 
space test for a complex object filling most of the inner 
grid will be run as a two-mesh problem, with monopole 
boundary conditions on the outermost surface. Complex tank 
test cases have been run with a truncated third grid to pro- 
vide the proper tank aspect ratio. 

It is recommended that, for all but the very simplest 
objects, the object definition input be developed piece by 
piece using NASCAP's OBJDEF-only mode. The main structure 
should first be plotted on quadrille paper and defined. 
NASCAP's printed and graphical output can then serve as an 
aid to adding various features and appendages. Complex ap- 
pendages (such as curved solar panels) can be defined 
separately and later added to the main object definition 
file. Finally, the defined object should be carefully checked 
before proceeding with the simulation. 

1.4.2 Flux De f inition; Operating Modes 

The calculation of the f.'.ux to the object dictates the 
mode of operation of the CHARGE portion of tJASCAP. The three 
flux-types are (1) the nonoenergetic electron gun source; -* 

(2) the Maxwellian probe approximation; and (3) reverse- 
trajectory sar.pl ing. 

For the monoenergstic electron source the user is re- 
quired to provide the current density as measured on a cali- 
bration plate. 1JASCAP automatically determines the number 
of particles to be tracked forward to reproduce the input 
calibration to reasonable accuracy. The user should verify 
the accuracy using the printed and graphical current density 
output. NASCAP automatically compensates for the presence of 
a constant magnetic field. 
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The Maxwellian probe approximation is the simplest and 
fastest operating mode. The environment is specified by the 
four parameters of ion and electron temperature and density. 
Since the angular distribution (isotropic) can be treated 
analytically, electron backscatter and low energy electron 
emission reduce to one-dimensional integrals. The most severe 
shortcoming of this approximation is the neglect of particle 
shadowing. 

The reverse-trajectory particle-pushing mode is the 
slowest operating mode for NASCAP. This is because good 
resolution requires a reasonable density of particles on the 
three-dimensional space of energy, polar angle, and azimuthal 
angle, resulting in a very large number of particles. Given 
the user-specified numbers of energies and angles, NASCAP 
chooses their values in an optimized manner. The distribu- 
tion sampled Can be either a Maxwellian or can be taken from 
experimental data. Data from selected ATS-5 observations 
are provided with NASCAP. 

1.4.3 Floating, Fixed, Biased and Coupled Conductors 

The object defined by NASCAP may be composed of up to 
seven conducting sections. By default, all conductors are 
floating and are coupled only through the fields in the empty 
space surrounding the object. The "Keyword” input provides a 
medium for altering these specifications. 

Typically, the main coupling between two conductors is 
the capacitance through a glue joint or dielectric spacer. 
NASCAP accepts this information through the "CIJ" keyword. 

Any conductor may be held at fixed potential relative 
to plasma ground or test tank ground. The potential solver 
will then change the charge stored on that conductor as neces- 
sary. Since the total charge on the object is not held 
constant, one should not use monopole boundary conditions 
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(see next section) . Holding a conductor at fixed potential 
improves the convergence rate of the potential solver. 

Conductors may be biased relative to conductor 1, which 
is considered to be spacecraft ground. This results in redis- 
tribution of charge among the conductors. Monopole boundary 
conditions may still be used. 

1.4.4 Use of the Potential Solver; Initial Potentials; 

Potential Scaling 

WASCAP's potential solver is discussed in detail in 
Chapter 7. As with any iterative scheme for solving differen- 
tial equations, the first few iterations consist of rapid local 
adjustment to minimize the residuals, followed by slow global 
drift toward the correct solution. By providing a good initial 
guess upon entering the potential solver, the second process 
becomes less important than the first, and thus the required 
number of iterations is minimized. The efficiency of the 
potential solver is decreased as dielectrics and vacuum gaps 
are made -thinner; it is increased when potentials on conductors 
are fixed. 

At the start of a simulation the potential may be set 
to zero everywhere (using IOUTER = 0) or to a monopole (Q/r) 
potential (using IOUTER = 2) . If a previous set of potentials 
is available for the same object, they may be scaled for the 
new run (using IOUTER = -n) . Since the potentials on the outer 
surface of the computational space are never changed by the 
potential solver, the initial guess also serves to define the 
boundary conditions. We recommend IOUTER = 0 for the test 
tank case, and IOUTER = 2 for the space cases. 

NASCAP provides the capability to scale the potential 
at the beginning of a simulation, and at each cycle prior to 
entering the potential solver. By default, the potential is 
scaled in proportion to the total charge ("SCALE") unless 
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conductor 1 is fixed ("NOSCALE"). The third option ("DSCALE") 
scales the potential according to 

(total charge) - (charge on fixed conductors) . 

This option is recommended for the test tank with a grounded 
conductor . 

We recommend that the number of potential iterations, 
MAXITR, be in the range 10 <_ MAXITR 50. The first ten 
iterations usually eliminate most of the error in the poten- 
tials. Beyond fifty iterations the algorithm becomes subject 
to roundoff error. Large numbers of iterations are usually 
required only in the early stage of a simulation. Here we 
suggest either a short time step, a potential restart, or 
the temporary fixing of conductor potentials. For an object 
with no differential charging, no potential iterations are 
required after the initial cycle, since the potential scales 
with the total charge. 

1.4.5 Choice of Time Step; LONGTIMESTEP Option 

The timescales involved in spacecraft charging range 
-5 —3 

from 10 - 10 seconds for formation of a photosheath about 

— 3 — 1 

a positive object, through 10 - 10 * seconds for an object 

to charge negatively in a magnetospheric substorm, to 1 - 
1000 seconds for differential charging. At the beginning 
of the simulation a fairly short time stop should be taken 
until a satisfactory potential solution is achieved. The 
time step can then be lengthened as the net current decreases. 
In a multicycle run, this can be done through the input varia- 
ble DELFAC . 

The explicit timestepping proceuure becomes unstable 

when 


dJ > 
dQ * 
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where J is the net current, Q is the total charge, and 6 is 
the time step. This situation is likely to happen when an 
object is charging positively (so that the derivative is 
large) or during differential charging (when a large 6 is 
desirable). For these cases the LONGTIMESTEP option should 
be invoked. This option causes the fluxes at the advance 
time step to be calculated using scaled potentials, and, if 
necessary, the current fluxes corrected using a first-order 
implicit scheme. The use of the LONGTIMESTEP option is not 
recommended in the presence of fixed cond’ ctor potentials, 
or when the "SHEATH" option is invoked. 

1.4.6 Treatment of Low Energy Electrons; The "SHEATH" Option 

In the default ("NOSHEATH") mode low energy emitted 
electrons (secondaries and photoelectrons) are treated simply 
as positive incident currents if the emitting surface is 
negative. For an electron-attractive emitting surface, this 
current is reduced by exp (-E*?/2) where E is the field at the 
emitting surface, is the outward normal with a length of one 
grid unit, and the characteristic emitted energy is 2 eV. 

If the "SHEATH" option is invoked, the electrons con- 
stitute g the currents described above are tracked, starting 
with 2 eV energy normal to the emitting surfaces. The space 
charge associated with these electrons is calculated for use 
by the pote itial solver. The associated currents are used 
to correct the fluxes to each surface cell, thus taking into 
account any photocurrent from one part of the satellite to 
another. 
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2. NASCAP RUNDECK AND INPUT SPECIFICATIONS 


This chapter contains specifications and examples for 
the NASCAP rundeck and three of the four user-supplied input 
files required by NASCAP. The specifications for the object 
definition input appear in Chapter 6. 

2.1 NASCAP FILES AND RUNDECK 

2.1.1 NASCAP Files 

NASCAP uses a substantial number of external files. 
Since the major portion of the rundeck is devoted to the 
assignment/ allocation/ and/or copying of these files, it is 
appropriate to discuss them prior to discussing the rundeck. 

NASCAP external files {Table 2.1) may be divided into 
three main categories: 

1. User-generated input files. These are the flux 
and object definition files, the keyword input 
file, and the runstream input. 

2. Scratch files. 

3. Code-generated files required for restart. These 
contain potentials, charge distributions, and 
object definition information. (If no restart 

is contemplated, they may be regarded as scratch 
files. The IAREA file may be regarded as a 
restart file if successive runs are done wi:h 
the same sun orientation.) 

EXEC 8 and NASCAP necessarily communicate file references by 
Logical Unit Number (LUN) or dataset reference number (DSRN) . 

3 

The values suggested in Table 2.1 are those to which we at S 
have become accustomed, and will be used in all examples. 
However, the LUN's are all input values and thus at the user's 
discretion. 
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2.1.2 The NASCAP Rundeck 

The HASCAP rundeck contains the EXEC 8 system commands 
required for execution of HASCAP. It performs (in order) the 
following functions: 

1. The @RUN card sets priority, charge number, time 
limit, etc. 

2. Necessary file assignments and copy operations are 
performed . 

3. A copy of the NASCAP absolute element is secured. 

4. The system is instructed to execute HASCAP. 

5. The runstream input file is supplied. 

6. A postmortem dump may be requested. 

7. Copy operations maybe performed to prepare for 
subsequent restart. 

In the examples which follow, we adopt the following 
conventions : 

1. All LUN's take the values suggested in rra bJ g f, ^- 1 - 

2. Catalogued files are assigned for exclusive use 
(@ASG,AX) . Conversely, 'all files so assigned are 
assumed to be catalogued. 

3. NASCAP communicates only with user files and 
scratch files. Restart files are copied into 
scratch files prior to execution and updated 
afterward. (This prevents their loss due to ab- 
normal termination.) Such catalogued files are 
named 10A., etc. 

4. The NASCAP absolute element is assumed to exist 
under che name NPROG.PGM. 
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5. The runstream, flux definition, object definition, 
and keyword input data are denoted [RUNSTREAM] , 
[FLUX], [OBJECT], [KEYWORD] respectively. Files 
containing these data are assumed named RUNSTREAM. , 
FLUX., OBJECT., and KEYWORD. 
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A self-contained runstream involving neither 
restart nor catalogued input file. 
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Example 2. A runstream appropriate to restart and using 
catalogued input files. 
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Example 3. The same as example 2, except 

1. The DeForest environmental data is placed on 
scratch file 9 for NASCAP use, and 

2. The IAREA file is regarded as a restart file. 
(See Section 2.1.1.) 
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2.2 


RUNSTREAM INPUT 


NASCAP runstream input, which may follow the @XQT 
statement in the run stream, contains grid size specifica- 
tions, logical unit numbers for all required files, restart 
(see Table 2.2) and plot flags, boundary condition and flux 
flags, potential solution limits and flags, time cycle 
limits and increment, sun direction vector and intensity, 
and, finally, flags and specifications for any optional plots 
desired. 

A list of all variables read from this data set ap- 
pears as Table 2.3. Unless otherwise noted, all variables 
are one per card, with integers in FORMAT (15) and real 
numbers in FORMAT (F10.0). The remainder of this section 
consists of the runstream input specifications (Table 2.3), 
and several runstream examples. 
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Table 2.3.' Runstreara Input Variables 


FORTRAN 

Variable 

Card No. Name Description 


GRID DEFINITION 


1 NX 

2 NY 

3 HZ 

4 NG 


Number of x-airection grid points. 
The only acceptable entry is 17. 

Number of y-direction grid points. 
The only acceptable entry is 17. 

Number of z-direction grid points. 
Hust be of the form 4n+l, with 

4 1 n 1 8 * 

Number of nested NX x NY x NZ grids. 
Potential solving is done in 1IG 
grids; particle tracking is done in 
at most two grids. 


5 


6 


7 


8 


9 


10 


LOGICAL UNIT NUMBERS (LUN’s) 


IP 


IAUN 


IR 


IY 


IU 


ISPARE 


LUN of external file on which grid 
potentials will be stored (and on 
which potentials from previous run 
lie if IPUEST =1). 

LUN of external file on which co- 
efficient product (Au) n values for 
grid points will be stored. 

LUN of external file on which con- 
jugate gradient array r values for 
grid points will be stored. 

LUN of external file on which 
the conjugate gradient array 
y values for grid points will 
be stored. 

LUN of external file on which the 
conjugate gradient array u values 
for grid points will be stored. 

LUN of the ‘'spare" external file 
which will rotate values with IP# 
IR, IY, and IU . 
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Table 2.3 

Card No. 
11 

12 

13 

14 

15 


16 

17 

18 


i 


. Continued 

FORTRAN 

Variable 

Name Description 

IROUS LUN of the external file on which 

the ROUS charge array will be 
stored. 


IPCOND 


IQCOND 


I OBJ 


IOBPLT 


I SAT 


LUN of the external file onto which 
the PCGND array will be dumped upon 
run completion for restart purposes. 

LUN of the external file onto which 
the QCOND array, ICYCS, and TIME 
will be dumped upon run completion 
for restart purposes. 

LUN of external file on which OBJDEF 
will write the PIICND, LIST, WPCOUD, 
and W arrays for use by OBJSPC. 

LUN of external file on which OBJDEF 
will write satellite block and cell 
definition data. 

IOBPLT is also a flag indicating 
whether or not satellite illustra- 
tion plots are desired. 

IOBPLT > 0 implies plot material 
composition satellite silhouettes, 
non-hidden line satellite building- 
block plots, and hidden-line cell- 
by-cell plots will be plotted. 

IOBPLT < 0 implies do not plot satel- 
lite illustration plots, and use unit 
number IPLOT I ABS (IOBPLT) . 

LUN of the object definition input 
file. 


ICNOW LUN of the external file on which 

the CHOW array (current surface 
charge) will be stored. 

IFLUX LUN of the flux definition input 

file. 
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Table 2.3. Continued 


FORTRAN 

Variable 


Card No. 

Name 



Description 

19 

IKEYWD 

LUN 

of 

the keyword input file. 

20 

IAREA 

LUN 

of 

external file on which 


shadowed surface cell areas are 
stored. 


If IAREA > 0, a new sun direction 
vector, SUN (3), is used to make a 
new shadowing calculation, and the 
resultant areas will be stored on 
unit IAREA. 

If IAREA < 0, the shadowed areas 
from a previous run will be used and 
are presumed to be stored already 
on unit IABS (IAREA). 

21 IPART LUN of external file onto which PUSH 

will write particle trajectory data 
for plotting by PARPLT. 

Also a flag indicating whether or 
not particle trajectory plots (for 
non-tank test cases) are desired. 

IPART > 0 implies write particle 
trajectory data for particles 
landing on cell IOCL(l) onto file 
IPART, and plot these trajectories 
for the first time cycle. 


RESTART, ODJDEF , PLOT, BOUNDARY CONDITION 
AND FLUX OPTIONS 


22 

ICREST 

Restart flag for restarting at a new 
time cycle. See Table 2.2. 

23 

IPREST 

Restart flag for potential solution. 
See Table 2.2. 
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Table 2.3 

Card No. 
24 


25 


26 


27 


23 


t 


I 


. Continued 


FORTRAN 

Variable 

Name 


IOBCAL 


ICON 


ITPART 


ITCUR 


IROUSP 


Description 

Flag for OBJDEF call. 

IOBCAL = 1 implies call OBJDEF and 
proceed normally. 

IOBCAL = 0 implies OBJDEF output is 
already available on files IOBJ and 
IABS (IOBPLT) . 

IOBCAL = -1 implies call OBJDEF (and 
SATPLT and KIDCEL if called for) and 
stop. 

Potential contour plot flag indicating 
whether or not and how often potential 
contour plots are desired. 

ICON <. 0 implies no plots. 

ICON > 0 implies plot potential con- 
tours every ICON time cycles, start- 
ing with cycle 1. 

Plot flag for tank test particle 
trajectory plots. 

ITPART = 0 xmplies no plots. 

ITPART r 0 implies plot particle 
trajectories. 

Plot flao for tank test current con- 
tour plots. 

ITCUR = 0 inplies no plots. 

ITCUR / 0 implies plot current con- 
tours . 

Plot flag for ROUS (space charge den- 
sity) contour plot. 

IROUSP = 0 implies no plots. 

IROUSP 5 i 0 implies plot final ROUS 
contours. 
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Table 2.3. Continued 


FORTRAN 

Variable 

Card No . Wane 


29 ICNVP 


30 IOUTER 


31 ITYPE 


DescriDtion 


Potential convergence printer plot 
flag indicating whether or not and 
how often printer plots of con- 
vergence data are desired. 

ICNVP = 0 implies no plots. 

ICNVP > 0 implies plot convergence 
data every ICNVP tine cycles starting 
with cycle 1. 

Boundary condition and potential 
initial guess option. 

IOUTER = 0 implies begin with zero 
potential everywhere. 

IOUTER = 1 implies begin with 1/r 
potentials on outer boundary of sys- 
tem; zeroes in all other locations. 

IOUTER = 2 implies begin with 1/r 
potentials everywhere. 

IOUTER = -N (negative integer) im- 
plies multiply previous potentials 
and charges by lo -1 ' if ICREST = 0 
and IPREST = 1. f Start at cycle 1 

with scaled version of previously 
calculated potentials. (See Table 2.2) 

Flux definition flag (must agree with 
first card of flux definition file) . 

ITYPE = 1 implies test tank case. 

ITYPE = 2 implies isotropic Maxwellian 
plasma approximation. 

ITYPE = 3 implies particle-pushing 
space case (experimental data or 
Maxwellian ambient plasma) . 


POTENTIAL CONVERGENCE PARAMETERS 

32 IALG Flag indicating which potential solu- 

tion algorithm to use (0 or 1). (See 
Section 7.1.) Normally IALG = 1. 
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Table 2.3 

Card No. 

33 

34 


35 


36 

37 

38 

39 

40 


B 


. Continued 


FORTRAN 

Variable 

Nane 


MAX IT R 


INTITR 


POTEPS 


Description 

Nunber of potential iterations de- 
sired. Zero is allowed; 50 is a 
suggested maximum. (Note possibility 
of potential restart.) 

Flag indicating that the conjugate 
gradient algorithm is to be re- 
initialized every INTITR iterations. 
INTITR is not normally used and 
should be set to a large nunber (e.g., 
999). 

Intended as an epsilon for determin- 
ing quality of potential convergence. 
Not currently used. 


TIME CYCLE INFORMATION 


NCYC 


Number of time cycles to be run. 


MCYC 


DELTA 

DELFAC 


POTENT call flag. 

MCYC = 0 implies POTENT will not be 
called to calculate new potentials; 
scaled potentials will be used. 

MCYC > 0 implies POTENT will be 
called every MCYC time cycles begin- 
ning with cycle 1. 

Length of time step in seconds for 
one charging cycle (time cycle) . 

Factor by which DELTA is multiplied 
at each cycle. 


SUN DEFINITION 

SUN (3) FORMAT (3F10.0). Vector indicating 

direction from center of inner grid 
to sun for shadowing calculation. 
SUN will be used to calculate new 
* shadowed areas only if IAREA > 0. 

SUN need not be normalized. 
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Table 2 . 3 

Card No. 

41 

42 


43a 


» 


- ! 


1 


. Continued 


FORTRAN 

Variable 

Name Description 

SFAC Solar intensity relative to 1 a.u. 

solar constant. 


REQUESTS FOR OPTIONAL PLOTS 


NCON, 
NDIR, 
NDIV (6) 


I 


READ 3020/ NCON, NDIR, (1JDIV(I),I = 1,6) 

3020 FORMAT (815) 

NCON is the number of additional 
potential contour plot cuts desired. 

NDIR is the number of additional 
viewing directions for which the 
3-D satellite illustration plots are 
desired. 

NDIV (I) is the number of extra divi- 
sion plots desired for the Ith view 
of the material composition silhouette 
plots (views 1 through 6 correspond 
to +x, -x, +v, -y, +z, and -z) . 

If this card contains all zeroes (or 
is blank) , no further cards will be 
read in. 

If NCON is non-zero, we read: 


ICNOPT (2,8) READ 3010, ( (ICNOPT (I , J) , I = 1,2),' 
J = 1, NCON) 

3010 FORMAT (1615) 

ICNOPT (1,J) = 1, 2, or 3 indicating 
the Jth cut is through fixed y, 
or z plane. 

ICNOPT (2, J) is an integer between 1 
and NX, NY, or NZ (depending on 
ICNOPT (1 , J) , indicating the fixed 
value of the cut. This value is in 
units of the next-to-inner grid. 

If NDIR is non-zero*, we read: 
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V I 

Table 2.3 

Card No. 
43b 


43c 


i 


l 


Continued 


FORTRAN 

Variable 

Name 

DIROPT (3,5) 


NDVAL 

(2,5,6) 


Description 

DO 100 I = 1, NDIR 

READ 1030, (DIROPT (J , I) , J = 1,3) 

100 CONTINUE 

1030 FORMAT (3F10.0) 

DIROPT (1,1), DIROPT (2,1), and 
DIROPT (3,1) form the direction-of- 
view vector from the satellite cen- 
ter to the viewer for the Ith 
optional view. DIROPT need not be 
normalized . 

For each non-zero entry of NDIV, we 
read one card: 

DO 200 IVIEVJ = 1,6 

NDIVI = NDIV(IVIEW) 

If (NDIVI. EQ.0) go to 200 

READ 3020, IVIEVJ2 , ( (NDVAL(I,J,IVIEW2) , 
I = 1,2) , J = 1, _ NDIVI) 

If (IVIEU2.NE. IVIEW) RETURN 0 

3020 FORMAT (15, 5X, 1015) 

200 CONTINUE 

For a given view IVIEW, a card is 
read in with a check variable IVIEW2 
(which must equal IVIEW) and NDVAL. 
NDVAL (1,J,IVIEW2) and NDVAL (2, J, IVIEW2) 
are the lower and upper Units of the 
Jth optional division for view IVIEW. 
The optional plot for this division 
will plot shaded surface areas for 
only those cells which lie within 
these limits (limits for the default 
plots are 1 to NX, 1 to NY, and 1 to 
NZ) . 
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4 (Runstream 

Input) . 

A 

particle-pushxn< 

run is to be 

re c tarted 

for three cycles. 

Is contained 

in one cubic 

grid. No 

plot: 


SOH'OIRECIIOm 


The object 


quested. The potentials are to be scaled at each 
cycle (no CG iterations) . The first cycle run will 
be 0.5 seconds, and the time step will increase by 
30 percent each cycle. The object is in the dark. 
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Example 5 

(Runstream Input) . 

A space 

case 

in the Maxwell probe 

approximation is to be 

run with 

two 

full grids, starting 

from time 

zero. The object has 

been defined previously. 

and 

a set 

of potential 

s exist;. 

The area calculation 


will be performed for the sun shining from the (111) 
grid direction. Object plots and potential contour 
plots will be made. One cycle will be run for 
5 x 10“5 seconds, including 25 Cj potential iterations. 
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Example 6 (Runstream Input) . A test tank case is to be run 
in three full grids. This is a new problem. The 
object will be defined and object plots made. Tank 
particle trajectory plots and current density contour 
plots are requested. One cycle will be run fcr 
5 x 10~5 seconds, with 25 potential iterations. 
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2.3 


FLUX DEFINITION INPUT 


The first card of the IFLUX file must contain the word 
'TYPE' in columns 1-4 and 1, 2, or 3 in column 10. The value 
must agree with card 31 of the runstream input. The types are: 

1 - Test tank, electron beam case 

2 - Isotropic Maxwellian space case 

3 - Particle-pushing space case 

Subsequent cards differ for each case,. 


Type 1 - Test tank with electron beam source 


Card 2 

Card 3 

Card 4 

Card 5 
Card 6 

Card 7 

Card 8 

Card 9 

Cards 10- 
45 


- IBMCAL (14) 

Normally 0; set to 1 for stop after 
initial beam calibration. 

- 2-location of beam calibration plane 
in grid units (F10.1) . 

- Beam size (X and Y) at the calibration 
plane (meters) (2F10.2). 

- Be=>m energy (eV) - (F10.1). 

- Initial electron velocity in code units 
(E10.2) . 

- Gun location in outer grid units (3F10.1). 

- Sample pJane (in grid units) (F10.1). 

- Total beam current (amperes) (F10.1). 

- Current density calibration values (11F4.1). 


All beam information is calibrated using the maximum 
value of the beam size from card 4. The current density array 
contains the current density profile as a function of (R,0) 
for 36 evenly spaced angles and 11 evenly spaced radial values 
Each card read contains the 11 radiai values at a single 
angular value. 
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x beam size (card 4) 


r 


i 


(i-1? 

10 


6 . 


l 


(i-1) 

36 


x 2 tt 


The central value is repeated on each card. The values are 
read in arbitrary units and the code will normalize the total 
beam current to the input value. Card 8, the sample plane 
location should be on or behind the target; it is used for 
terminating particle trajectories. (See Figures 2.1, 2.2.) 


Type 2 - Maxwellian probe approximation 

All cards are FORMAT (F10.0, A6) . Accepted literals 
are MKS and CGS for plasma density, EV, JOULES, KEV, and 
KELVIN for plasma temperature. The first density/energy re- 
fers to electrons and the last to protons. The literal END 
(optional) is also accepted, and terminates the flux definition. 

Type 3 - Flux determined by revf rse trajectory particle tracking 

Card 2 contains one of the two literals (A6) MAXWEL or 
DEFOR. Following 'MAXWEL' input is as in Type 2, (q.v.), 
except an END card is required. Following ' DEFCR’ , cards are 
required giving IUNIT, RKOUR, IDAY (I5/F1 0 . 0/15 ) . The data 
contained in the ELT NASCAP.DEFOR must have been copied onto 
file IUNIT. 

Table 2.4 lists the input RHOUR, IDAY for which data 
from ATS-5 is in the DEFOR file. The particle spectra for 
these times are described in Appendix A of Reference 1. Any 
other data may be used if it is prepared in the following 
fashion. The data must consist of 63 card images with the 
following format: 
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» Nt • « 


X 


4 = 


Figure 


4 = 0 



2.2. Specification of electron gun characteristics for 
the test tank case. Electron fluxes at the sample 
plane for an incharged environment are to be 
specified at she solid points. 
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Table 2.4. Hourj and Days for ATS-5 Flux Dat^ Contained 
in File Element NASCAP .DEFOR. 


Date 

RHOUR 

IDAY 

3/14/71 

5.999 

73 


7.000 

73 


7.501 

73 


8.001 

73 


8.502 

73 


9.998 

73 

12/3/70 

1.997 

337 


3.999 

337 


6.001 

337 


7.002 

337 


7.997 

337 


12.001 

337 

3/18/70 

8.998 

77 


9.999 

77 


10.818 

77 


11.205 

77 


12.001 

77 


17.006 

77 
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5 FORMAT (12, IX, 4F12.2, 2F8.3, 14, F4.3) 
The Read list is 


READ (IUNIT, 5) K, (FDEF ( J , K) , J=l,5), RHOUR, IDAY, AVMIN 

K = card number (1-63) and energy bin number. 

2 

FDEF (1,K) = parallel electron flux in eV/(cm sec sr eV) . 

2 

FDEF (2 , K) = parallel proton flux in eV/(cn sec sr eV) . 

FDEF (3, K) = perpendicular electron flux in eV/ 

2 

(cm sec sr eV) . 

FDEF (4, K) = perpendicular proton flux in eV/ 

2 

(cm sec sr eV) . 

FDEF (5, K) = lower bound on energy bin (eV) . 

RHOUR = checked against requested hour. 

IDAY = checked against requested day. 

AVMIN = not used. 


Missing flux data is indicated with a 0. FSPACE assumes that 
there is no flux above 50 keV and that data on cards past 28 
are for energies above 1000 eV. The flux for a velocity V 
is assumed to be 


ly-Ysi 

f<V) ~ |V| |V g | f parallel 




f 


perpendicular 


where V s is a direction vector, whose FORTRAN name is STV 
(see below) . 

Subsequent cards contain 


NSPEC 

(15) 

Number of species 

NENG 

(15) 

Number of energy es 

NTHET 

(15) 

Number of polar angles 

NPHI 

(15) 

Number of azimuthal 
angles 

NSTP 

(15) 

Maximum number of steps 
per particle 
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RMASS(l), CHARGE (1) 


(2F10.16) 

Mass (kg) , charge (C) 

RMASS (NSPEC) , CHARGE (NSPEC) 

VCODE (F10.6) Initial code velocity 

STV(3) (3F10.0) Vector characterizing 

flux anisotropy 

The present code requires NSPEC = 2, with species 1. 
being electrons and species 2 protons. The final card, STV, 
may be omitted, in which case the flux is assumed isotropic. 
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Example 7 (Flux Definition) . Defines a neutral 10 keV 

plasma for use in the Maxwell probe operating mode. 
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Example 9 (Flux Definition) . Takes flux from DeForest 'en- 
vironmental data for hour 9.998 of day 73. The data 
is read from LUN 9. 
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Example 10 (Flux Definition). Defines a non-neutral, non- 
equilibrium Maxwellian plasma for use by the reverse 
trajectory particle pushing mode. 


40 



ORKIINJAI. l> 

of poor ^a^ity 


ItPE 1 

0 

21.0 

.30 

20U0.0 

• sotoo 

».o 

2 I .0 

a. tut-; 


10. 

io. or 

• u 

9 

• 0 

a 

3 

7 

• 6 

6 

6 

4 

• 7 

3 

2 

i 

.5 

0.5 

0 

1 1 . 

10.0 9 

• 0 

9 

• c 

8 

C 

7 

.-9 

6 

S 

9h 5 

2 

9 

i 

.9 

C.S 

0 

12. 

1 0 . 0 9 


9 

• 0 

fl 

0 

6 

• 9 

b 

b 

4 

4 

2 

0 

i 

. 9 

0.5 

0 

13. 

10. 0 v 

* 6 

ti 

• 9 

7 

6 

6 

« d 

b 

6 

4 

4 

2 

6 

i 

.3 

0.5 

0 

11. 

10.0 9 

• 6 

tt 

• b 

7 

5 

6 

• J 

b 

3 

3 

V 

2 

3 

0 

.9 

0.5 

0 

1%. 

10.0 9 

« 7 

6 

• 7 

7 

9 

6 

.0 

1 

S 

3 

u 

1 

7 

0 

.5 

0.0 

c 

1 6 . 

10.0 9 

• 0 

b 

• 6 

7 

5 

5 

.7 

3 

9 

2 

7 

1 

0 

0 

.5 

C.O 

0 

17. 

1C.0 9 

• y 

6 

.7 

7 

3 

5 

*b 

3 

4 


b 

0 

b 

0 

.5 

0.0 

0 

18. 

1C .0 1 u 

• 0 

tt 

• 9 

7 

9 

S 

« b 

3 

3 


7 

u 

3 

0 

.5 

0.0 

0 

19. 

10. C 9 

• i 

b 

.7 

7 

2 

S 

• 4 

3 

2 


S 

C 

1 

0 

.5 

0.0 

3 

20. 

10.0 9 

• 0 

e 

.5 

6 

a 

4 

. V 

2 

7 


1 

0 

b 

0 

• C 

C.O 

0 

21 . 

1L.0 9 

• 8 

6 

• 4 

6 

J 

4 

• b 

2 

b 


0 

c 

b 

0 

.0 

0.0 

0 

22. 

10.0 V 

• e 

6 

• S 

6 

7 

4 

.7 

2 

7 


I 

u 

5 

0 

.0 

0.0 

0 

23. 

IC.O v 

.9 

0 

• S 

6 

V 

S 

• C 

3 

t 


4 

0 

1 

0 

.5 

0.0 

3 

29. 

IC.O 9 

. V 

8 

• 6 

7 

1 

S 

.2 

3 

b 


7 

u 

6 

c 

.5 

J.O 

0 

25. 

1 u . 0 1 0 

• 0 

6 

.7 

7 

9 

5 

• 6 

3 

9 

2 

6 

1 

1 

0 

. 1 

0.5 

0 

26. 

1 J.01O 

• t 

9 

• 0 

7 

7 

6 

• 3 

4 

8 

3 

L 

l 

0 

0 

.9 

0.5 

0 

27. 

10.010 

• 0 

9 


8 

3 

7 

• b 

b 

b 

4 

a 

2 

b 

0 

. 7 

0.5 

0 

28. 

1C. Clio 

• u 

9 

• 6 

9 

0 

6 

• & 

7 

b 

s 

b 

J 

b 

I 

.5 

0.5 

0 

29. 

1C. 01c 

* c.* 

V 

.9 

9 

c 

b 

• ; 

6 

6 

s 

tt 

4 

0 

2 

. 1 

0.5 

3 

30. 

I0.CIC 

• c 

9 

• 6 

9 

0 

8 

. v 

6 

V 

s 

U 

J 

a 

1 

. V 

0.5 

0 

31. 

IQ.tJU 

• L' 

9 

• 4 

9 

0 

8 

• / 

b 

6 

b 

4 

3 

b 

l 

. 6 

0.5 

0 

12. 

10.01 J 

• 0 

y 

• 4 

9 

0 

7 

• b 

6 

6 

H 

7 

3 

C 

1 

. 9 

0.5 

u 

J3« 

1 e . 0 1 e 

• 0 

9 

• 4 

t> 

7 

7 

■ 6 

& 

& 

3 

6 

l 

2 

1 

.0 

0.5 

3 

39, 

10. Ole 

• 0 

9 

• 3 

8 

5 

6 

• 6 

S 

2 

3 

U 

l 

6 

0 

.5 

vi . 0 

0 

35. 

1C. Old 

*C 

V 

• 3 

o 

1 

& 

• 2 

4 

4 

2 

b 

1 

2 

0 

. 2 

0.5 

0 

36. 

10.01 0 

• 0 

V 

• 4 

a 

2 

6 

• 3 

4 

3 

2 

4 

i 

1 

0 

. 1 

0.5 

u 

37. 

13.01. 

• c 

9 

• b 

b 

3 

b 

• 4 

4 

* 

2 

4 

i 

1 

0 

. 1 

e . 5 

0 

38. 

1 0 . u 1 e 

• c 

9 

. 4 

8 

1 

6 


4 

3 

2 

5 

1 

2 

0 

. 2 

0.5 

z 

39. 

10.010 

.0 

9 

.3 

b 

l 

6 

• 2 

4 

4 

2 

t 

1 

b 

u 

.5 

3.0 

0 

90. 

1 J.CI 3 

• c 

9 

.2 

8 

2 

6 

• 3 

4 

o 

3 

H 

2 

0 

0 

.a 

0.5 

3 

9 t . 

I 1. . 0 t w 

• c« 

9 

• 1 

6 

1 

6 

• 0 

b 

& 

4 

U 

2 

4 

I 

. i 

0.5 

0 

‘■2. 

1 e . 0 1 e 

• c 

V 

• 1 

8 

2 

7 

. 3 

b 

3 

4 

b 

2 

7 

I 

.5 

0.5 

3 

93. 

1..V/I. 

• -/ 

9 

.2 

o 

3 

7 

. 4 

6 

e 

b 

0 

J 

1 

1 

• 9 

j . s 

C 

99 • 

1 . • 0 1 . 

• c 

V 

• l 

8 

H 

7 

6 ^ 

6 

b 

4 

7 

J 

4 

1 

.5 

u.5 

a 

95. 

1 ..Ole 

• 1 

9 

• C 

o 

9 

7 

• 6 

6 

b 

4 

6 

3 

J 

4 

.5 

0.0 

0 


Example 11 (Flux Definition). Defines flux in the tank test 
mode with a 2 keV beam. 
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2.4 


KEYWORD INPUT 


Additional input to NASCAP is contained in file IKEYWD. 
This file must be present and have proper attributes, though 
it may contain no card images. The format is a 'keyword* 
beginning in column 1, of which the first six characters are 
significant, and data (if any) beginning in column 21. With 
the exceptions noted below, cards may appear in any order. 



Keyword 

Data 

Data Format 

1. 

COMMENT 

none 


2. 

END 

none 


3. 

SURFACE CELL 

ICELL 

(110) 

4. 

LONGTIMESTEP 

none 


5. 

SHEATH 

none 


6. 

NOSHEATH 

none 


7. 

SCALE 

none 


8. 

NOSCALE 

none 


9. 

DSCALE 

none 


10. 

TANKSIZE 

IL, IR 

(215) 

11. 

BFIELD 

B ( 3 ) 

(3F10.0) 

12. 

DIPOLE 

P<3) , R (3 ) 

(6F10.0) 

13. 

CIJ 

I, J, C 

(215, F10.0) 

14. 

BIAS 

I, VBIAS 

(110, F10.0) 

15. 

FIXP 

I, VFIX 

(110, F10 . 0 ) 
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1 . 

2 . 

3 . 

4 . 

5 . -6 
7.-8 

9 . 

10 . 

11 . 

12 . 

13 . 

14 . 

15 . 


\ 

r 


EXPLANATIONS 


Ignore card. 

Ignore subsequent cards. 

Print detailed information for surface cell numbered 
ICELL. If particle plots are made, they are for the 
lowest numbered such cell. 

Evoke the LONGTIMESTEP algorithm. This is useful for 
an object whose monopole potential is fully developed 
but is being differentially charged. 

. Treat space charge and fluxes due to emitted low energy 
electrons to first order. The default is NOSHEATH. 

. The potential may be scaled from the previous result 
according to the total charge prior to the conjugate 
gradient calculation. The default is SCALE if con- 
ductor 1 is floating and NOSCALE if conductor 1 is 
fixed . 

Potential to be scaled according to [total charge - 
(sum of charge on fixed potential conductors)]. 

The Z-dimension of the outer grid may be truncated from 
the left and/or right. (The outer grid must include 
the next inner grid.) 

2 

Constant magnetic field in Webers/m . 

•*■7 -*• 

Magnetic dipole with moment P(A - n>') at location R 

(grid units) . 

Capacitance C (farads) between conductors I and J. 

This is to be used for conductors which are glued 
together or are separated by a thin insulating layer; 
capacitance between disjoint conductors is handled 
implicitly. 

Conductor I to be biased at VBIAS (volts) relative 
to conductor 1, which is considered spacecraft ground. 
Biasing must be done in order, beginning with 
conductor 2 . 

Potential of conductor I to be fixed at ''T’lX. Space- 
craft ground may be fixed. 
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SURE ACE CELL 


55 


2. 

Cl J 

1 

2 

J97.L- 

3. 

SURE ACE CELL 


216 


•a. 

SOW ACt CELL 


29b 


5. 

SoWALC CLLL 


312 


6 • 

SoWACE CELL 


9 36 


7. 

noSCAlE 




a. 

r ixp 


1 

-575. 

*. 

SOW ACE cell 


9 


10. 

SuK» ACC CELL 


92 


1 1. 

SuWACE CELL 


95 


12. 

SOWACE CELL 


7b 


13. 

Sow ACt CCLL 


31 


19. 

SuWACE CELL 


10J 


IS. 

SOW Ut cell 


15 


1 A . 

S u w ace cell 


9 


17. 

TAUAS I 2t 

09 

19 


18. 

E*<0 




19. 

3f 1 EL U 

0. 


.000052 

20. 

E-.D 




21 . 

a I AS 


2 

70. 

22. 

E-.O 




23. 

CIO 

2 

3 

.00005 

29. 

CIO 

1 

3 

.009 

25 • 

BIAS 


2 

10. 

26. 

Fit P 


0 

33. 

27. 

E»«0 





-.000019 


Example 12 (Keyword Input) . Several surface cells are chosen 
for output. A capac'tance of 347 pf is placed between 
conductors 1 and 2. Potential scaling is suppressed. 
Conductor 1 is fixed at -57S volts. The outermost 
grid is truncated at subscripts 2=4 and 2 = 14. 
Remaining cards are ignored (although in correct format) . 
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t . 

SJHTaCL CELL 

2. 

0SC«Lt 

3. 

r i xp 

9 . 

SUMTACC CELL 

s. 

Suh» *CE CELL 

4. 

L0N<* f t ME STEP 

7. 

I.JSriLATH 

8 . 

TANAS 1 ZE 

9. 

8f 1 LLD 

10. 

END 


S3 

1 03. 

9 
92 


oa 27 

o. .003CS2 -.U00019 


Example 13 (Keyword Input) . Three surface cells are chosen 
for output. The DSCALE potential scaling opt: on is 
invoked. Conductor 1 is grounded. The LONGT1 MESTEP 
option is invoked and the SHEATH option is suppressec . 
The outer grid is truncated at subscripts Z = 3 and 
Z = 27. A constant magnetic field is present. 


4b , 




3 . OUTPUT 


3.1 GENERAL CONSIDERATIONS 

NASCAP provides a great deal of output, both printed 
and graphical, to aid the user in following the progress of 
the simulation. To a large extent NASCAP output is self- 
explanatory. In this chapter we provide further description 
to help the user understand and anticipate NASCAP output. 

3.2 PRINTED OUTPUT 

NASCAP printed output includes: 

1. Runstream input echo (self-explanatory) 

2. Object definition output (Section 6.3) 

3. Shadowing calculation output (Section 3.2.1) 

4. Flux definition output (self-explanatory) 

5. Keyword input echo and summary (Section 3.2.2) 

6. Time cycle output 

a. , CHARGE segment output (Section 3.2.3) 

b. POTENT segment output (Section 3.2.4) 

7. Final values of potential array (Section 3.2.5) 

In addition, self-explanatory messages are printed at various 
points to help the user monitor the progress of the cede and 
determine whether the desired actions are being taken. 

3.2.1 Shadowing Calculation P rintout 

Subroutine HIDCEL prints out a variety of diagnostic 
information when HIDCEL is called for a shadowing calcula- 
tion. (When called for satellite illustration, HIDCEL 
prints only the final number of shielded A1 polygons, NA1.) 
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HIDCEL first prints a list of the initial unshielded 
A1 polygons: For each polygon, the following data is 

printed : 

IA1: Index of this A1 polygon. 

IA1T: Index of surface cell which generated this 

A1 polygon. 

NVM1: Number of vertices of this A1 r . Ton. 

AREA2(IALT): Unshielded area of this A xygon. 

Next, the number of initial A1 polygons, NA1, the 
number of initial A2 polygons, NA2, and a print flag, IP, 
are printed. 

Next, the set of final A1 polygons and their final 
(shielded) areas are printed. The surface cell index IA1T 
of each polygon's parent cell is also printed. If a poly- 
gon occurs which contributes to the area of a surface cell 
which an earlier polygon (or polygons) has contributed to, 
a message to this effect is printed. 

A list of the final fractional areas (= original 
unshielded area/final shielded area) for each surface cell 
is printed, along with the original unshielded and final 
shielded areas themselves. 

Finally, the set of shielded A1 polygons in 2-D is 
explicitly listed by printing the complete vertex set for 
each polygon (including the duplicate first vertex at the 
end of the .list) . 

3.2.2 Keyword Input Echo and Summary 

Each card read by the keyword input routine, INKEYW, 
is echoed. In addition, many of the keywords elicit an im- 
mediate response verifying the requested action. The mes- 
sage "DECODE ERROR" results from an improper data format, 
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which results in program termination. An unrecognized 
keyword also results in termination. 

Following the input echo, a summary is presented. 

The surface cells specified for flux I/O are listed. The 
code units of charge and capacitance are given. (V = Q/4rR, 
where V is in volts and Q and R are in code units.) The 
matrix of capacitances between conductors is printed. The 
values of the LONGTIMESTEP, SHEATH, and potential scaling 
options are given. Truncation of the outer grid is noted. 
The constant magnetic field and any magnetic dipoles present 
are listed. 

3.2.3 CHARGE Segment Output 

The CHARG3 segment printed output is as follows: 

1. Following a page feed, the cycle number, time 
total charge, and potentials and charges on 
conductors are printed. 

2. For the test tank case, some of the internal 
variables are printed, followed by the array 
CURENT (3, NX, NY) , the incident current at the 
sample plane. 

3. Detailed information is printed for each sur- 
face cell specified in the keyword input. 

This includes the interna] code defining the 
surface cell, its location, outward normal , 
and composition, the mean potential on the 
surface and the mean normal el-ectric field 

in the thin dielectric, the various incident 
and emitted flux components, and the net flux. 

The potentials, fields, and fluxes are those at 
the beginning time of the time step. 

4. The net charging current to the object is 
printed. 
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5. If the LONGTIMESTEP option has been invoked, 
items 3 and 4 are repeated for scaled poten- 
tials at the advance time. 

6. If the SHEATH option has been invoked, the 
corrected fluxes are printed for each of the 
selected surface cells. 

7. The updated time is printed, followed by the 
CNOW array (the charge on each surface cell) . 

3.2.4 Potential Solution Printout 

QSUM, the sum of all charges in the system is printed 
out initially, as are the PCOND and QCOND arrays. NC, the 
number of conductors in this problem is printed. Next, 

MBIAS, the number of biased conductors (1 if none are 
biased) and the VBIAS array containing the amounts of volt- 
age bias are printed. 

Next, several variables to be used by OBJSPC are 
printed: SWPCIJ, SSCIJ, SWPSSC, arrays CIJSUM, SWPCND, 

and ROUSCN. 

Finally, we enter the potential solution proper, and 
the initial value of RDOTR is printed. The potential itera- 
tion loop prints out a repeated series of convergence re- 
lated data including: iteration number, QCOND (1), SWPCND (1), 

AUNC(l), PCOND(l), RCOND(l), UCOND(l), the PCOND and QCOND 
arrays, UDOTAU (u n « (Au) n ), ALPHA, RDOTR and RDOTRS (r n+1 -r n+1 
and r n *r n ), BETA, C and CS <c n+1 and c n ) , RTDOT (r»r, the 
monotonic residual), and, finally, UDOTU (u *u ). 

Upon completion of tne potential solution and call 
to QCONCP, the final PCOND and QCOND arrays are again 
printed. 
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If convergence printer plot data was called for by 
ICNVP > 0, printer plots follow: the first plots the log 

of RDOTR versus iteration number; the second plots either 
(algorithm 0) the log of UDOTAU versus iteration number, 
or (algorithm 1) the monotonic residual RTDOT against 
iteration number; the third plots the potential on con- 
ductor one versus iteration number. 

3.2.5 Potential Array Printout 

Final potentials are the last output for any run in 
which POTENT has been called. The successive blocks of 
numbers may be regarded as X-Y cuts through a potential map. 
Successive grids of data are printed, beginning with the 
innermost grid. For the innermost grid, zeroes indicate 
the object interior. For outer grids, the next inner grid 
appears as a block of zeroes. 


3.3 


PLOT OUTPUT FROM NASCAP 


NASCAP generates a wide variety of diagnostic plots. 

Any of these may be requested or suppressed via user-defined 
input flags in the main input file. 

3.3.1 Satellite Illustration Plots 

Several types of plots illustrating the satellite are 
created when IOBPLT > 0. These plots serve not only as 
report quality output, but also as diagnostic output for 
the programmer himself. Defining a complex satellite in 
three-dimensions, complete with a variety of surface materials 
is a task prone to errors, especially for those unfamiliar 
with the code. The satellite illustration plots are an im- 
portant early-on indication of whether or not the satellite 
as imagined has been prox^erly communicated to the code. 

The first of these plots are the material composition 
silhouettes. From each of the six possible on-axis directions, 
a view of the satellite is plotted showing each surface cell 
and the material it is comprised of via up to fifteen dif- 
ferent area shadings. In a complex satellite involving some 
inward-facing surfaces, surface cells hidden by another sur- 
face will not be shown by the six -default views. However, 
the user may request additional plots for any or all of the 
six viewing directions specifying a given subset (or subsets) 
of the possible x, y, or z values such that only those surface 
cells falling within the specified range will be plotte.d. 

In this way a set of silhouettes may be generated which il- 
lustrates the material composition of all surfaces of tne 
satellite. 

Two types of three-dimensional perspective plots aid 
in picturing the satellite: (1) non-hidden-line perspective 

plots of each of the satellite "blocks” as defined by the user 
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in input file ISAT; and (2) hidden-line perspective plots of 
the satellite showing all surface cell boundaries. The first 
plot echoes back the user's definition of the satellite form, 
and the second illustrates the three-dimensional shape and code 
resolution of the satellite in a clear fashion. Each of these 
plots is plotted for each of three different default viewing 
directions, viz. (.5,. 8, .5), (-.5,. 8,. 5), and (.2, -.5, .5). The 
user may also specify up to five additional viewing directions. 

An additional hidden-line cell-by-cell perspective plot 
will be plotted as a diagnostic feature when a surface cell 
shadowing calculation is requested {I AREA > 0) . In this case 
the direction-of-view is defined by the SUN input vector. 

3.3.2 Potential Contour Plots 

Potential contour plots will be generated when ICON > 0. 
Contours may be plotted along any cut across the inner t*.o 
grids; however, the default action is to plot contours along 
the x-y, x-z, and y-z planes through the center of the grid. 

Up to eight additional cuts may be requested by user input. 

The projection of the satellite perimeter is also plotted 
with the contours. (Calling for potential contours for a one- 
grid problem will result in abnormal termination.) 

3.3.3 ROUS Contour Plo ts 

Contours of the ROUS charge array as of the final cycle 
will be plotted when IROUS > 0. Contour plot cuts are the 
same as the default cuts for potential contours: the x-y, 

x-z , and y-z planes passing vhrough the center of the grid. 
However, since ROUS is defined only for the inner grid, only 
one grid is plotted. 

The proiection of the satellite perimeter is again 
plotted with the contours. 
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3.3.4 Particle Trajectory Plots (Non-Tank Test Case) 

When IPART > 0, particle trajectories across the inner 
two grids will be plotted for the first time cycle. Three 
trajectory projections (x-y, x-z, and y-z planes) will be 
plotted for each of the two particle species (electrons and 
protons). Only trajectory data for those particles which 
hit the lowest numbered cell specified in the KEYWORD input 
(or cell 1 if none have been specified) are plotted. 

The projection of the satellite perimeter is also 
plotted with the trajectories. 

3.3.5 Particle Trajectory Plot (Tank Test) 

When ITPART j* 0, particle trajectories projected onto 
the x-z plane will be plotted. 

The projection of the satellite perimeter will also 
be plotted. 

3.3.6 Current Contour Plot (Tank Test) 

When ITCUR ^ 0, current contours along the x-y plane 
through the center of the grid wixl be plotted. 

The projection of the satellite perimeter will also 
be plotted. 
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4 . SAMPLE RUN 


As a crude illustration of how NASCAP may be used, a 
short sample run is presented. This run begins a simulation 
of a teflon coated sphere with nominal diameter 0.9 meters 
in a neutral 14 keV plasma of density 1 cm \ The flux 
definition, keyword and ob 3 ect definition inputs are shown 
in Figures 4.1 through 4.3. 

The printed run output appears in pages 57 to 91. 

Prior to the first cycle appear the runstream input echo, 
the object definition output (Section 6.3), the shadowing 
calculation output (Section 3.2.1), the flux definition echo, 
and the keyword input echo and summary. 

Cycle 1 begins on page 68. The negligible initial 
charges and potentials result from the 1/4 rrr boundary con- 
dition and initial potential guess. The problem is slightly 
asymmetric because it is not possible to locate a sphere of 
diameter 3 mesh units exactly at the center of the computa- 
tional space. For this case, the net flux to the sphere is 
*v»7 percent of the incident electron flux. 

The potential calculation may be traced through the 
printout on pages 69 through 74, and is summarised by the 
plots on pages 75 through 77. The non-monotonic residual 
(RDOTR) and the decreasing reconstructed residual (RTDCT, 
although the title is UJOTAU) illustrate the rapid conver- 
gence during the early iterations, followed by slower con- 
vergence in the later cycles. Though the solution is con- 
tinuing to improve after the twenty- five requested itera- 
tions, the plot of the conductor potential (third plot) 
indicates that further iterations are not worthwhile. 

The second cycle (page 78) begins at 0.001 seconds 
with a potential of -11 volts. This is sufficient to 
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Figure 4.1 

. Object definition input for sample 

run. 



TTF*l 2 

i.oa cos 

i.o cos 

J *1 « JO KE ¥ 

M.O KEN 


Figure 4.2. Flux definition input for sample run. 


1. SURFACE OIL l* 

2. SURF * CE LLLL ■» 

3. END 


Figure 4.3. Keyword input for sample run. 


55 



X 


> 



1 



! 

i 

I 


noticeably reduce the net charging current. At the end of 
the simulation (0.0024 seconds) the sphere is at a poten- 
tial of 25 volts. 

Finally, the potential array at the end of the 
simulation is printed. The graphical output from this run 
(pages 92 to 95) shows the quasisphere as viewed from the 
sun direction and three sets of potential contours cal- 
culated after the cycle 1 conjugate gradient iteration pro- 
cess. On each potential contour plot the minimum and maxi- 
mum potentials in volts (ZMIN, ZMAX) in the plane, and the 
interval between contours (AZ) are indicated. 
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5. MAIN, MISCELLANEOUS AND UTILITY ROUTINES 


In this section are described the main program, 
NASCAP, the sub-main program, TRILIN, their subsidiary 
routines, and various utility routines used in the 
NASCAP code. 

5.1 MAJOR ROUTINES (NASCAP, TRILIN) 

NASCAP 

Purpose ; NASCAP is the main routine. All other subroutines 
are called either directly or indirectly by NASCAP. 

Internal Organization ; (See flowchart. Figure 5.1.) 

NASCAP calls INOPT to read the runstream input file de- 
fining the grid, logical unit numbers and execution options. 

If the input options call for any of the available 
plots, a flag is set to specify opening and closing of the 
plct output file at the appropriate time. 

IOBCAL is checked; If IOBCAL is non-zero, OBJDEF is 
called to read sate) lite definition input and define the 
satellite. If IOBCAL is zero, the satellite has been de- 
fined by a previous run and OBJDEF is not called. In this 
case, OEJDEF output data required at this point is read from 
file IPLOT = IABS (IOBPLT) , and this file is rewound. 

If any plots were called for, a plot file is opened. 
If satellite illustration plots were requested, SATPLT is 
called to generate all three types of plots. 

It a shadowed area calculation is requested for this 
run (IAREA > 0) , HIDCEL is called to calculate the areas 
and generate the corresponding hidden-line plot. 

If IOBCAL equals -1, this run is solely for object 
definition and plotting; solution of the problem is skipped. 
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figure 5.1. NASCAP flow chart 
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Otherwise, subroutine TRILIN is called to cycle through 
the appropriate number of time steps, charging up the 
satellite, push. ng particles, and solving potentials. 

Upon return from IRILIN, if a plot file has been 
opened, it is closed and execution terminates. 


TRILIN 

Internal Organization : (See flowchart. Figure 5.2.) After 

calling FLXDEF and INDATA for initialization purposes, serves 
as primary program for NASCAP main loop: alternate calls to 

CHARGE and POTENT. Calls necessary routines for potential 
scaling, potential convergence plots and potential contour 
plots. After end of main loop, calls routines for restart 
dump, particle trajectory plots and charge density contour 
plots. 

Called 3y : NASCAP 

5.2 I/O ROUTINES 
APRT 

Purpose : APRT prints NGP grids' worth of data for the 

NX x NY x N2 x NG array stored on file IF. 

Calling Sequence : SUBROUTINE APRT (BUF1, IF, NGP, NX, NY, 

NZ , NG) 

Called By : POTENT (most calls axe currently commented out) , 

INDATA 

Local Variables : IF: Unit number of file containing data 

to be printed. 

NGP: Number of grids' worth of this data 

to be printed. 
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FLXDEF 

Purpose : Reads and preprocesses flux definition parameters. 

Calling Sequence : SUBROUTINE FLXDEF (ITYPE1, NX, NY, NZ, NG) 

Called By : TRILIN 

GETNC 

Purpose : Extract NC, the highest conductor index in problem, 

from JSURF (NSURF ) . 

Calling Sequence : SUBROUTINE GETNC (NC) 

Called By : INDATA 

INDATA 

Purpose : Performs data manipulations associated with 

initialization/restart. 

Calling Sequence : SUBROUTINE INDATA (NX, NY, NZ , NG, NC) 

Called By : TRILIN 

INKEYW 

Purpose : Reads "keyword" input file and sets up associated 

options. 

Calling Sequence : SUBROUTINE INKEYW (IKEY) 

IKEY - LUN of keyword input file. 

Called Bv: INDATA 
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LORDER 


Purpose : Puts in numerical order LIST(l) entries starting 

at LIST (2) . If W(l) f 0.0, array W parallels LIST(2)... 

LIST (LIST (1}+1) and will be reordered similarly. (Identi- 
cal routine LORDEQ is called by OBJDEF. 

Calling Sequence : SUBROUTINE LORDER (LIST, W) 

Called By : INKEYW 

INOPT 

Pu rpose : INOPT reads the runstream input file (corresponding 

to card input) which defines the grid limits, logical unit 
numbers, restart, plot and execution options, potential iter- 
ation limits, the time cycle limits and timestep. 

INOPT prints out all of these data at the start of 
NASCAP printed output. 

Calling Sequence : SUBROUTINE INOPT (NX, NY, NZ, NG, NC, 

MAX IDF, PEPSDF ) 

Called By : NASCAP 


I 
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MOVDAT 

Purpose ; MOVDAT is the general data movement subroutine 
used throughout the code to clock transfer large quanti- 
ties of data between external storage files and internal 
code storage arrays. 


Calling Sequence : SUBROUTINE MOVDAT (N, DAT, IUNIT, IN) 


Called By : 


(All routines 

doing block-transfer 

I/O) 

CHARGE 

PUPDAT 

QSUMER 

HIDCEL 

RUPDAT 

SETALL 

OBJDEF 

UDOTUC 


POTENT 

YUPDAT 


TRILIN 

UUPDAT 


ROUSP 

UUPDAU 


OBJSPC 

QCONCP 


URSETO 

SETOP 


YINIT 

PMODQ 



Local Variables; 


N: Number of words to be transferred. 

DAT: Core storage array data is to be 

moved into or out of. 


IUNIT: Logical unit number of file to be 

read or written onto. 

IN: IN = 1 implies read from file into 

DAT; IN = 0 implies write from DAT 
onto file. 
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REWIND 

Purpose : REWIND rewinds file IUNIT by executing- the block- 

transfer (NTRAN) rewind command. For block-transfer files, 
the regular FORTRAN rewind cannot be used. 

Calling Sequence ; SUBROUTINE REWIND (IUNIT) 

Called By : (All routines doing block-transfer I/O) 


CHARGE 

RUPDAT 

HIDCEL 

UDOTUC 

OBJDEF 

YUPDAT 

POTENT 

UUPDAT 

TRILIN 

UUPDAU 

ROUSP 

CCONCP 

OBJSPC 

SETOP 

URSETO 

PMODQ 

YIN IT 

QSUMER 

PUPDAT 

SETALL 


Local Variables : IUNIT: Logical unit to be rewound. 

RSDUMP (See INDATA) 

Purpose : Prepares files for possiole future restart. 

Calling Sequence :-- SUBROUTINE INDATA (NX, NY, NZ , NG, NC) 

ENTRY RSDUMP (NX, NY, NZ , NG, NC) 

Called Bv: TRILIN 


SUMQD 

Purpose : Used for "DSCALE" option. Calculates (and prints) 

OSUMD - QSUM - f £ d CCOND(I fixed ) 
conductors 

Calling Sequence : SUBROUTINE SUMQD (QSUMD , QSUM, QCOND, NC) 

Cal led By : TRILIN 
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5.3 POTENTIAL SCALING AND RELATED ROUTINES 
PMODQ 

Purpose : Scales potentials in accordance with total charge 

on object. 

Calling Sequence : SUBROUTINE PMODQ (BUF1, IP, ISPARE, PCOND, 

NX, NY, NS, NG, QSUMO, QSUM) 

Called By : TRILIN 

QSUMER 

Purpose : Calculates total charge contained within inner grid. 

Calling Sequence : SUBROUTINE QSUMER (IROUS, ROUS, QCOND , NX, 

NY, NZ, QSUM) 

Called By : TRILIN, INDATA 

RESETQ 

Purpose : Scales all potentials and charges by 10**IOUTER. 

Called ior ICREST = 0, IPREST = 1, IOUTER < 0. 

Calling Sequence : SUBROUTINE RESETQ (NX, NY, NZ , NG, Q) 

Called L y: INDATA 


SETALL 

Purpose : Called when IPREST = ICREST = 0, IOUTER = 2 to 

set all potentials to Q/4~r. 

Calling Sequence : SUBROUTINE SETALL (NX, NY, NZ, NG, Q) 

Called Bv: .TNDATA 
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SET IP 

Purpose : For IOUTEF. = 0 or 1, SETIP writes the initial 

guess potentials for the innermost grid (zeroes, in these 
cases) onto file IP. 

Calling Sequence : SUBROUTINE SETIP (BUF1, IP, NX, NY, NZ, NG) 

Called By : INDATA 

SETOP 

Purpose : After the innermost grid data for the initial 

guess potentials have been placed on file IP by SETIP, SETOP 
is called to fill in the data for all remaining grids. 

Grids 2 through NG-1 are zero-filled, and the outer- 
most grid is either zero-filled (IOUTER = 0), or filled witn 
the l/4*r data generated by SETPOP for its outer edges 
(IOUTER = 1) . 

Calling Sequence : SUBROUTINE SETOP (BUF1, BUF2 , IP, NX, NY, 

NZ, NG) 

Called By ; TRILIN 
INDATA 

SETOP is called only when IOUTER = 0 or 1. 
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SETPOP 

Purpose : For IOUTER =1, SETPOP defines the initial guess 

potentials on the outer edges of the outer grid according 
to Coulomb's law: 

P = Q/4rr 

Potentials not on the outer edge of the outer grid 
are set to zero for this case. 

Calling Sequence : SUBROUTINE SETPOP (P, QSUM, NX, NY, NZ) 

Called By : TRILIN 

Local Variables : P: Storage area BUF1, representing outer 

grid initial guess potentials. 


5.4 CNVPLT AND PRINTER-PLOTTER 
CNVPLT 

Purpose : Plots potential convergence parameters using 

printer-plotter routines, i.e., ALOG10(R(I)) versus I, 
ALOGIO (U (I) ) versus I, P(I) versus I. 

Cal Lng Sequence : SUBROUTINE CNVPIT (R, U, P, NI) 

R(NI) - An array of positive numbers. 

U (NI) - An array of positive numbers. 

P (NI ) - An array assumed to contain 

potentials of spacecraft ground. 

NI - Dimension of arrays R, U, P. 

Called By : TRILIN 
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DOPLOT 


Purpose : Generates printer plot from data in COMMON blocks 

/PLOT/, /I0031/. (Part of printer plotter package.) 


Calling Sequence : SUBROUTINE DOPLOT (I5YM) 

ISYM(l) - array of symbols to be used in 
plot. 


Called By : CNVPLT 


INPLOT 


Purpose : Initialises printer plot. (Printer plotter 

package. ) 


Calling Sequence : 


SUBROUTINE INPLOT (XMIN, XMAX, YMIN, YMAX, 
IBANDS) 


XMIN ) 
XMAX ( 
YMIN ( 
YMAX ' 


axis limits of plot 


IBANDS - length of plot in units of ten 
lines 


C alled By : CNVPLT 


LINPLT 

Purpo se: Printer plotter routine. Plots point (X,Y) with 

symbol NUMSYM and straight line to previous point with sym- 
bol LINSYH. 

Calling Sequence : SUBROUTINE LINPLT (X, Y, NUMSYM, LINSYM) 

Called B y: CNVPLT 
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OLDLIN (See LINPLT) 

Purpose : Printer plotter package. This entry not currently 

used by NASCAP. 

Calling Sequence : SUBROUTINE LINPLT (X, Y, NUMSYM, LINSYM) 

Entry OLDLIN (X, Y, NUMSYM, LINSYM) 


POINT 

Purpose : Pringer plotter package. Plots point (X,Y) with 

symbol NUMSYM. 

Calling Sequence : SUBROUTINE POINT (X, Y, NUMSYM) 

Called By : CNVPLT 

REPLOT (See INPLOT) 

Purpose : Printer plotter package. This entry not used in 

NASCAP . 

Calling Sequence : Entry REPLOT (XMIN, XMAX, YMIN, YMAX, 

IBANDS) 
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6. CODE DESCRIPTION - OBJECT DEFINITION SECTION 

6.1 INTRODUCTION AND GENERAL CONSIDERATIONS 

6.1.1 Function 

The OBJDEF section of the NASCAP code serves to trans- 
late relatively simple user specifications of geometrical 
shapes, surface composition and material properties into the 
far more complex information required for efficient operation 
of the POTENT, CHARGE, HIDCEL and SATPLT sections. The user 
input to OBJDEF is read from file ISAT, which is described in 
Section 6.2. In addition to printed output (described in 
Section 6.3), OBJDEF generates files IOBJ (for use by POTENT) 
and IOBPLT (for use by CHARGE, SATPLT and HIDCEL). The contents 
of these files is given in Section 6.4. 

6.1.2 Volume Cells 

The object to be simulated is defined on the inner 
cubic mesh. The "lower left-hand" corner of the inner mesh 
is at mesh point (1,1,1), while the "upper right-hand" corner 
is at (NX,NY,NZ). The object should be restricted to the 
space 

2 <_ x < NX - 1 
2 < y < NY - 1 
2 < z _< NZ - 1, 

i.e., it should not contact the edge of the inner grid. For 
convenience, coordinates are specified relative to an "OFFSET", 
which is initially set to the center of the mesh, and may be 
specified by the user using the "OFFSET" command (see below) . 

An "OFFSET" of (0,0,0) will allcw object definition in terms 
of "absolute" coordinates. A volume cell is referenced by 
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its lowest indexed corner. The physical distance correspond- 
ing to a grid unit is specified using the "MESH SIZE" command. 

The object (as well c.s the space external to the 
object) is built of the four building blocks shown in Fig- 
ure 6.1: the cube, the right prism or "WEDGE", the tetra- 

hedron, and the complement of a tetrahedron. Of course, the 
user is not required to define an object cell by cell, but in 
terms of larger geometrical components (Figures 6.2-6.13): 

1. The rectangular parallelepiped ("RECTAN”) , 
defined by a corner coordinate and three edge 
lengths. (Object i on Figures 6.2-6.13.) 

2. The right prism ("WEDGE"), defined by a corner 
coordinate, two lengths and an orientation. 

(Object 5 on Figures 6.2-6.13.) 

3. The tetrahedron ("TETRAH") , defined by a corner 
coordinate, a length and an orientation. 

(Object 3 on Figures 6.2-6.13.) 

4. The "FIL111" object, a special type of wedge 
described in Section 6.2.7. (Teflon-coated 
portion of object 6 on Figures 6.2-6.13.) 

5. The octagonal right cylinder ("OCTAGO"), a 
complex object formed of two " RECTAN "s and four 
"WEDGE "s , defined by the ends of its axis, its 
width and irs side. (Object 4 on Fiaures 6.2- 
6.13.) 

6. The quasisphere ( "QSPHER" ) , a multiply complex 
object formed of three "OCTAGO" s and eight 
"TETRAH "s , defined by its center, diameter and 
side. (Object 2 on Figures 6.2-6.13.) 

Detailed input specifications for each of the above objects 
are given in Section 6.2. 
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Figure 6.1. Four shapes of volume cells considered by the 
NASCAP code: (a) empty cube; (b) wedge-shaped 

cell with 110 surface; (c) tetrahedron with 111 
surface; (d) truncated cube with 111 surface. 
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Figure 6.2. Illustration of NASCAF object definition. Non- 
perspective view showing surface cell compositions. 
(1) rectangular parallelepiped; (2) quasi-sphere; 
(3) tetrahedron; (4) right octagonal cylinder; 

(5) right triangular prism; (6) compound object 
consisting of aluminum rectangular parallelepiped, 
aluminum right triangular prism, and teflon <111 > 
wedge . 
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Figure 6.3. Another view of Figure 6.2, 
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Figure 6.5. Another view of Figure 6.2. 
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6.1.3 Surface Cells 

The four types of volume cells give rise to five dis- 
tinct types of surface elements, or surface cells. A sur- 
face cell is identified by {1} the volume cell within which 
it is located, or out of which its normal points, and (2) 
the direction of its outward normal. The latter is conven- 
iently expressed in crystallographic notation: (100), (101), 

etc., and the divisions of these directions into different 
symmetry classes provide a primary categorization of surface 
cell types. The surface cell types are: 

(100): a. simple square, size lxl (mesh units); 

b. right isosceles triangle, which must be 
further identified by the location of 
its right-angle corner; 

(110) : rectangle 1 x /2 (mesh units) associated with 

wedge cell; 

(111) : equilateral triangle of side /l (mesh units) ; 

there are two types of such cells, depending 
on whether the associated volume cell is ia) 

1/6 full and 5/6 empty, or (b) 5/6 full and 
1/6 empty. 

In addition to its type and location, a surface cell is 
labeled as to material composition and as to the conductor it 
overlays. NASCAP's internal surface cell code format is shown 
in Figure 6.14. 

A maximum of 1024 surface cells is allowed by NASCAP. 
Definition of superfluous surface cells is automatically by- 
passed. Surface cells that become superfluous after their 
definition may be eliminated by the "COMPRE" command, initi- 
ating a call to subroutine CMPRSS. (CMPRSS is automatically 
called following the 'ENDSAT" command or an EOF on unit ISAT.) 
In case of a multiply - defined surface cell , the latest defini- 
tion supersedes the earlier. 
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5 43 210 987654 321098 765432 1C9876 5 43210 

U 1 1 ! II II II I I I I 

HGF E D C BA 


Field Bits 
A 4-0 

B 5 


C 11-6 

D 17-12 

E 23-18 

F 29-24 

G 32-30 

H 34-33 


Material index 

Set for right-triangular 100 sur- 
faces and for 111 surfaces whose 
enclosing volume cell is mostly 
empty. 

Direction of surface normal (in 
crystallographic notation) 

Z-coordinate 

Y-coordinate 

X-coordinate 

Conductor index 

Orientation code for right- 
triangular 100 surfaces. 


Figure 6.14. Surface element list (JSURF) entry format. 
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6.1.4 Special Cells, Surface Points, Special Points, 

Multiple Conductors 

A volume cell is labeled as a special cell if a poten- 
tial interpolation formula other than the usual trilinear 
one must be used; i.e., if it is partially filled or if a 
right triangle surface cell points into it. The element 
table, LTBL (NX,NY,NZ), (see Figure 6.15), is set up by the 
OBJDEF section for use elsewhere to identify special cells. 

A point is a surface point if it is a vertex of a sur- 
face cell. A surface point is a multiconductor point if it 
is a vertex of surface cells overlaying different conductors. 

Conductors are numbered 1 < IC < 7 using the "CONDUC" 
command (q.v.). Conductor 1 is usually considered spacecraft 
ground. (See Keyword Input, Section 2.4.) 

A vertex of a special cell which is not a surface point 
is a special point . 


6.1.5 Restrictions on Object Geometry 

1. No sjrface point or special point may lie on the 
mesh boundary. It follows that all surfaces must 

be at least one unit removed from the mesh boundary, 
and a right triangle surface cell must be at least 
two units from the mesh boundary toward which it 
points . 

2. A maximum of 1021 surface cells is allowed. 

3. No volume cell can have two or more right-triangle 
surface cells po '.nting into it. 

4. A partially filled volume cell cannot be super- 
seded by another partially filled cell except 
through use of the "DELETE" command (q.v.). 
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Code for Element Table [LTBL (NX,NY,NZ) ] . 


54321 C987654321 0 9 8 765 432109876 5 43210 

I 1 | I I 1 1 I 

E DC B A 


Field Bits 

A 4-0 

B 14-G 

C 18 

D 19 

E 30-21 


Cell-type code (see Appendix D) 

Orientation code (see Appendix D) 

Set if cell is completely filled 
(interior) 

Set for ar. empty special cell 

Index used to reference EKOJ array 
to determine low energy electron 
currents 


OAiEHTATION CODE 

3x3 bits. Each group of 3 contains 1, 2 or 3 in the 
lower 2 bits, with the high bit set for negative. 


Code 


takes 


m. m_ m_ 

(-) 1 i 1 , (-) i 2 , (-) 


r m. m_ m^ ■, 

(r ,r ,r ) to (-) 1 r. , (-) r. , (-) r. 

A J l X 1 X 2 3j 


e.g., the following codes take a point to (x,y,z): 


Octal Code 

123 

365 

532 

176 

567 

617 


Dec. Code 

1,2,3 
3, -2,-1 
-1,3,2 
1,-3, -2 
-1,-2, -3 
1,-3 


P oint 

(x,y,z) 

(-z,-y,x) 

(-x,z ,y) 

(x,z,-y) 

(-x,-y,-z) 

(y,-x,-z) 


Figure 6.15. Element table codes and orientation codes. 



6.2 


INPUT SPECIFICATIONS 


6.2.1 General Considerations 

User input to the object definition section of NASCAP 
is read from file ISAT primarily by subroutine INPUT and 
secondarily by subroutines RECTAN, OCTGON, WEDGE, -FIL111, 
QSPHER and TETRAH. The action taken by INPUT in response 
to a record is determined by a literal in columns 1-6. 
Recognized literals are indicated in Table 6.1. If the 
literal is an object name, one of the above subroutines will 
be called and will expect further input (see below) . Con- 
trol will be returned to INPUT upon encountering an "ENDOBJ" 
or "ENDSAT" card. If the literal is recognized to be a 
command, the corresponding data (if any) will appear on the 
same card. An unrecognized literal is assumed to be a 
material name and must be followed by three cards containing 
material parameters (FORMAT 8F10.0) (see Section 6.5). 

As an example, the input used to create Figures 6.2- 
6.13 is shown in Figure 6.16. 


6.2.2 Comment ("COMMEN") 

Card is ignored. 

6.2.3 Compress ("COMP RE”) 

Causes elimination of redundant surface cells. 

6.2.4 Conductor Index ("CONDUC") 

Causes subsequently defined surfaces to overlay the 
conductor whose index appears in column 15. The 
maximum conductor index is 7. The minimum and 
default is 1. 

6.2.5 Deletion ('TELETE") 

Causes the volume cell entered in columns 11-25 
of card (relative to current offset) [Format (315)] 
to be declared empty and all associated surface cells 
deleted. 
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Table 6.1. 

Literals Recognized by 
Commands 

Subroutine INPUT 

Literal <a) 

Data 

Data Format ^ 

COMMENT 

(none) 


COMPRESS 

(none) 


CONDUCTOR 

ICON 

(15) 

DELETE 

IJK 

(315) 

ENDOBJ 

(none) 

" 

END SAT 

(none) 


MESH SIZE 

XMESH 

(F10. 0) 

OFFSET 

IOFF 

(315) 



Objects 

' 


Literal (a) 

Subroutine Called 



FI1.111 

FILlll 



OCTAGON 

OCTGON 



QSPHERE 

QSPHER 



RECTANGULAR PARALLELEPIPED 

RECTAN 



TETRAHEDRON 

TETRAH 

- 

1 

WEDGE 

WEDGE 



* a ^Only first six characters 

significant . 



^ Begins in column 11. 
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ORIGINAL PAG l IS 
OF POOR QUALJi y 


1. 

aluminum 



2. 

1 . 

0.001 

-1 . 

3. 

2b0. 

1.73 

1.36 

b. 

17. 

18. 

19. 

S. 

TtFLON 



6 • 

2.000 

0.000127 

1 .t- 

7. 

2.2 

16.7 

1 .b 

a. 

17. 

ta. 

19. 

9. 

KAPTON 



10. 

3.503 

0.000127 

1 .£- 

1 1 . 

1 >b4 

9.3 

1 .b 

12. 

17. 

ia. 

19. 

13. 

MAGNES 1 UH 



lb. 

l • 

0.U01 

-1. 

IS. 

I. 7b 

2b. 32 

1 .36 

16. 

17. 

18. 

19. 

17. 

S I 02 



ta. 

b.000 

0.000127 

1.1- 

19. 

363. 

1 .63 

».b 

20. 

17. 

18. 

IV. 

21. 

CONOUC TOR 

1 


22. 

OFFSET 

2 

2 C2 

23. 

MESH SIZE 

• 025b 


2b. 

RECTANGULAR PARALLELEPIPED 

2S. 

corner 

0 

5 J 

26. 

DELTAS 

2 

} ! 

27. 

SURFACE 

♦ X 

ALUh I NUH 

48 • 

SuM» ACE 

.Y 

alum i N un 

29. 

SURFACE 

♦ 2 

AtUKlNUh 

3C. 

SURFACE 

-X 

A AH T ON 

31 . 

SURF ACE 

- r 

kahton 

32. 

surface 

-l 

KAPTON 

33. 

ENOOBJ 



3b. 

CONDUCTOR 

2 


3S. 

OFFSET 

12 

12 23 

36 . 

QSPnERE 



37. 

CiNTEK 

-3 

-2 0 

38. 

DIAMETER 

b 


39. 

SIDE 

2 


bO. 

MATERIAL 

teflon 


b 1 . 

ENOOBJ 



b2. 

C uNUUC T OR 

3 


b3. 

OFFSET 

12 

*4 2N 

bb. 

TlTkamEDRUN 


b 5 . 

CJHNLn 

-1 

-1 -I 

b 6 • 

r ace 

TEFLON 

1 

b 7 • 

Length 

2 


b a . 

S JRf ACE 

-X 

ALUM | NUH 

b 9 . 

SuKf ACE 

- y 

alum l num 

SO. 

5UHF ACE 

-2 

acum i Nun 

b 1 . 

E^OObJ 



52. 

COfti/UCTOri 

N 



13. 

•*0. 

40. 


lb 


10 . 


/ 0. 

43. 


lb 


70. 

40. 

14. 

<* 0 . 

20 . 


I b 


10. 


7U. 

20 . 


0.97 .3 

b.E-OS lb. 
21* 44. 

3. 0.3 

2. E-Ob 1 b . 


.24.0 . 

IS. 

23. 


■I. 

IS. 


21 . 

2. I 

.00002 

41 •' 

0-94 

•00004 

21 

2 .** 

.00002 

21 .* 


44. 

0. IS 

lb. 

22 . 

.2b 

lb. 

22 . 

O.b 

lb. 

22 . 


23. 


-I . 


2 . 


IS. 

23. 

2S0. 

IS. 

23. 


Object 1 

Rectangular parallepiped 


Obj act 2 

Quasi-sphere 


Object 3 

Tetrahedron 


Figure 6.16. Sample OBJDEF input. 


1 .3 
16. 
2b. 

0. 

16. 

2b. 

| 0 . 

16. 

2b. 

0. 

16. 

2b. 

. I 2 
*16* 
2M* 
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ORMPML PVjE lb 
OK POOR REALITY 


53. 

0* FSt I 


J 

1 2 22 

59. 

OC T AbUN 




56. 

AXIS 


CO 

t 

o 

o 

56. 

Mi 01H 


3 


57. 

s toe 


1 


5a. 

SURFACE 

• 


ALUMINUM 

59. 

SURFACE 

c 


5102 

60. 

E«0Ub J 




6 1 « 

CONDUCTOR 


5 


62. 

UFF5ET 


9 

6 12 

6 2. 

redue 




69. 

COnNEK 


-l 

-l -l 

65. 

r ACL 

TCFLOM 

0 

66. 

LLNbTH 


J 

2 2 

6 7. 

SURFACE 

-X 


ALUMINUM 

6 S . 

SURFACE 

- Y 


MAbNES 1 UM 

69. 

Surface 

-z 


ALUMINUM 

70. 

surface 

♦ X 


KAPTON 

71. 

C •JOOb J 




72. 

CuNOUCTOR 


6 


7 2. 

0FFSE7 


12 

12 2 

79. 

RtCTA.9bU|.AK 

parallelepiped 

75. 

C J9f.tR 


-1 

-1 0 

76. 

DELTAS 


9 

9 1 

77, 

SURF ACE 

♦ X 


ALUMINUM 

78, 

SURF ACE 

♦ Y 


aluminum 

79, 

surface 

♦ z 


ALUM 1NUM 

8u. 

Surf ACE 

•X 


ALUM 1 NUM 

81. 

Surface 

• r 


aluminum 

82. 

SURF AlC 

-z 


alum i nCim 

62. 

Euuuu J 




en. 

At jE 




65. 

C0-..9ER 


2 

2 1 

86. 

FACE 

ALUM1 N 

-1 

87. 

LlN.-TH 


2 

3 3 

88. 

SURF ACE 

♦ X 


ALUMINUM 

89. 

Surf ACE 

♦ Y 


aluminum 

9 Cl. 

SurF ACE 

-z 


aluc inum 

9 t. 

SjRF ACE 

♦ z 


aluminum 

92. 

CN0C8J 




92. 

offse r 


0 

0 0 

99. 

F IL 1 II 




96. 

C2RRER 


12 

15 J 

96. 

FACt 

TEFLON 

-l 

97. 

E .DOcJJ 




98. 

E.FDSA r 





i 


) 

) 


i 


i 


oo oo 


Object 4 

Right octagonal cy 


I Object 5 

Right triangular prism 


Object 6a 

Rectangular parallelepiped 


o 


Object 6b 

Right triangular prism 


12 3 ) Object 6c 

j <111 ;> wedge 


Figure 6.16. Sample OBJDDF input (continued). 
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6.2.6 End definition of geometrical component ("ENDOBJ" ) 
Control is returned to subroutine INPUT. 

6.2.7 End object definition ("ENDSAT") 

Control is returned to subroutine OBJDEF. 

6.2.8 Mesh Size ("MESH S") 

The mesh spacing (in meters) is taken to be the 
floating point number entered in columns 11-20. 

The default is 0.1 m. 

6.2.9 Origin Shift ("OFFSET”) 

Causes subsequent coordinate definitions to be taken 
relative to the absolute coordinates entered in 
columns 11-25 of card in Format (315). The default 
offset is the center of the mesh. 

6.2.10 <111> Wedge ("FILlll”) 

The <111> wedge is used for smoothing stepped surfaces. 
Its use is illustrated in object 6 of Figures 6. 2-6.1*. 

1. READ (11,21) CCODE , IC 

21 FORMAT (A6 , 4X , 315, 5X, 315) 

CCODE must contain the lateral "CORNER". 

IC (dimension (3,2)) contains the end coordinates 
of the right-angle-ccrner line. This line must be 
in a (110) symmetry direction. 

2. READ (11,402) CCODE, SMAT, ID 

Identical to description for "WEDGE", except 
(1) ID must be a (111) direction normal to 
line IC, and (2) all exposed surfaces will be 
of material SMAT. 
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6.2.11 Octagonal Right Cylinder ("OCTAGO") (see Figure 6.17). 

1. READ (11,21) CCODE , IC, ID 

21 FORMAT (A6, 4X, 315, 5X, 315) 

CCODE must contain the literal 'AXIS'. IC and ID 
are integer arrays containing the lower and upper 
indexed ends of the cylinder axis which must 
parallel a coordinate axis. If WIDTH is odd, the 
fraction one-half will be added in each direction 
normal to the axis. 

2.-3. READ (11,31) CCODE, ITMP 
31 FORMAT (A6 , 4X, 15) 

This read statement is executed twice. CCODE takes 
on the values 'WIDTH' and 'SIDE'. Integer variables 
of the same names are set to ITMP. It is required 
that WIDTH >_ 3. If SIDE >_ WIDTH, a value for SIDE 
will be provided. Otherwise, SIDE and WIDTH must 
be both even or both odd. 

4. READ (11,241) CCODE, ANS , SMAT 
241 FORMAT (A6, 4X , 2A1 , bX, A6) 

This input is similar to that for 'RECTAN' (6.2.13) 
except (1) only the first character of ANS is 
significant, the values '+', '-', and 'C' denoting, 
respectively, octagonal surfaces at the upper and 
lower ends of the axis and the circumferential sur- 
face, and (2) surfaces not explicitly defined are 
defaulted to the first material defined. 

6.2.12 Twenty-six Faceted Object ("QSPHER") 

1. READ (11,21) CCODE, CTR 

21 FORMAT (A6 , 4X, 315) 

CCODE must contain the literal "CENTER". 

CTR contains the coordinates of the Q-sphere 
center. For an odd diameter Q-sphere the true 
center will be (1/2, 1/2, 1/2) relative to the 
specified coordinates. 
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2. READ (11,241) CCODE, SMAT 
241 FORMAT (A6 , 4X, A6) 

CCODE may contain the literals "DIAMET" , "SIDE", 
"MATERI " , or "ENDOBJ" . "DIAMET" and "SIDE" are 
fully equivalent to "WIDTH" and "SIDE" for "OCTAGO" 
(q* v .); "SMAT" is decoded with format 15; the 
defaults are 3 and 1, respectivel-y . "MATERI" 
indicates the Q-sphere is covered with the material 
denoted by SMAT. 

6.2.13 Rectangular Parallelepiped ("RECTAN") 

1. READ (11,21) CCODE, IC 
21 FORMAT (A6 , 4X, 315) 

CCODE must contain the literal 'CORNER'. IC is an 
integer array containing the x, y and z coordinates 
of the lowest-indexed corner. 

2. READ (11,21) CCODE, ID 

CCODE must contain the literal 'DELTAS'. ID is an 
integer array containing the edge lengths parallel 
to the x, y and z axes. 

3. READ (11,241) CCODE, ANS, SMAT 
241 FORMAT (A6 , 4X, 2A1 , 8X, A6) 

CCODE may be 'SURFAC', 'COMMEN' , 'ENDOBJ' or 'ENDSAT'. 
In the case of a 'SURFAC' card, ANS (an array of 
length 2) contains the outward surface normal code 
in the form '+X', '-Y', etc. SMAT contains a 
material name which must have been previously 
defined. This read statement continues to be 
executed until an 'ENDOBJ' or 'ENDSAT' card is 
encountered, returning control to INPUT. 
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6.2.14 Tetrahedron ("TETRAH") 

1. READ (11,21) CCODE , IC 
21 FORMAT (A6, 4X, 315) 

CCODE must contain the literal "CORNER". 

IC contains the coordinate of the right-angle 
corner. 

2. READ (11,402) CCODE, SMAT, ID 
(See description for "WEDGE".) 

3. READ (11,21) CCODE, LENGTH 

(Same as for "WEDGE" except only the first 
value is significant.) 

4. READ (11,241) CCODE, ANS , SMAT 
(See description for "RECTAN" .) 

6.2.15 Right Triangular Right Cylinder ("WEDGE") 

1. READ (11,21) CCODE, IC 
21 FORMAT (A6, 4X, 315) 

CCODE must contain the literal 'CORNER.' IC is an 
integer array containing the x, y and z coordinates 
of the lower-indexed right-angle corner. 

2. READ (11,402) CCODE, SMAT, ID 
402 FORMAT (A6 , 4X, A6 , 4X, 315) 

CCODE must contain the literal 'FACE'. SMAT 
contains a material name which must have been 
previously defined. ID is an integer array 
containing the outward surface normal of the 
slanted face in crystallographic notation. 

3. READ (11,21) CCODE, LENGTH 

CCODE must contain the literal ’LENGTH’. LENGTH 
is an array containing the edge lengths parallel 
to the coordinate axes. 

4. READ (11,241) CCODE, ANS, SMAT 
(See description for ’RECTAN’.) 
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6.3 PRINTED OUTPUT 

Printed output from the object definition section con- 
sists of (1) modified input echo, (2) surface cell list, 

(3) potential coefficients for surface points when determined 
to be in error, (4) Element Table (LTBL) entries for special 
cells. Figure 6.15, (5) number of surface points, multi- 
conductor points and special points, (6) potential coefficients 
for special points when determined ;o be in error, and (7) 
tables of material properties. Also printed are statements of 
the form "CALLING subroutinenane" for diagnostic purposes. The 
output from defining the objects of Figures 6.2-6.13 is shown 
in Section 6.3.8. 

6.3.1 Modified Input Echo 

This portion of the output is self-explanatory. 

6.3.2 Surface Cell Li st 

A table of surface cells is printed by subroutine 
DIAGNO, including 

1. The cell number (array index). 

2. Tne internal code (octal) for the cell 
characteristics. 

3. The low-index corner (IX,IY,IZ) of the volume 
cell within which the surface lies, or out of 
which it poincs. 

4. The outward normal of the surface. 

5. The material name corresponding to the surface 
composition . 

6.3.3 Erroneous Potential Coefficients (Surface Points) 

Printed when a constant potential is not a local solution. 
(This out'put should not occur.) 
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6.3.4 Element Table (LTBL) Entries for Special Cells 
Printed for diagnostic purposes (see Figure 6.15). 

6.3.5 Numbers of Surface, Multiconductor and Special Points 

This output is self-explanatory. The maximum allowed 
values are 1024, 100 and 1024 , respectively. 

6.3.6 Erroneous Potential Coefficients (Special Points) 

(See Section 6.3.3). 

6.3.7 Material Properties 

The input values of the various material properties 
and their conversion to code units (see Section 6.5). 

6.3.8 Printed Output from Figure 6.16 
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6.4 


OUTPUT TO FILES IOBPLT, IOBJ 


6.4.1 File I PLOT [=ABS (IOBPLT) ] 

IPLOT is a FORTRAN unformatted file. The contents of 
successive records are listed below. Quantities in parentheses 
are dimensions. See glossary for definitions. 

1. XMESH 

2. NMAT, MCOD£(15), MATPR<20,15) 

3. NSURF , JSURF (1024 ) 

4. NBLK, NPBLK (64), IPER (3 , 32 , 64) , IFAC (80) 

5. NSPTS, NPTS 

6. SWPCND (7 ) 

7. LTBL (17,17,33) 

6.4.2 File IOBJ 

IOBJ is an NT RAN file intended for use by routine 
OBJSPC. It contains 

1. Potential coefficients linking nulticonductor 
points to conductors (P21CND) 

2. List of surface points 

3. Potential coefficients linking surface points to 
conductors 

4. Space potential coefficients for surface points 

5. List of special points 

6. Space potential coefficients for special points 
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6-5 MATERIAL PROPERTIES 

t 




< 




6.5.1 Definition 


The present version of NASCAP makes use of thirteen 
material properties; storage is provided for twenty. The 
twenty values are read on three cards following the material 
name using FORMAT (F10.0). The input values are converted 
to code values by subroutine MATPRO and its subsidiary rou- 
tines . 


Property 1: 
Property 2 : 

Property 3 : 

Property 4 : 
Property 5 : 

Property 6 : 


Relative dielectric constant (dimensionless) . 

Thickness of dielectric film or vacuum gap 
(meters) . 

c-lectrical conductivity (mho/m) . The value 
- 1 . indicates a vacuum gap over a conducting 
surface . 

Atomic number (dimensionless) . 

Maximum secondary electron yield for electron 
impact at normal incidence (dimensionless) . 

Primary electron energy to produce maximum 
yield at normal incidence (keV) . 


Properties 7-10; Range for incident electrons. Either : 

P r 

Range = P ? E 8 + ? g E ' 0 

where the 'range is m angstroms and for the 
energy m keV, 

or 

P 7 = -1. to indicate use of Feldman* s^ 
empirical range formula 

P g = density (g/cm 8 ) 

P^q = mean atomic weight (dimensionless) 

Property 11: Secondary electron yield for normally incident 

1 keV protons. 
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Property 12: Proton energy to produce maximum secondary 

electron yield (keV) . 

Property 13: Photoelectron yield for normally incident 

sunlight (A/m2) . 


6.5.2 Sample Values 

- 

There is . 

a great deal of unce 

rtainty in all property 

values concerned 

with particle range 

, low energy electron 

emission or conductivity of insulators. With that warning, 

we present below 

the input values we 

have tentatively adopted 

for clean aluminum, clean magnesium. 

t 

teflon, kapton and SiC^. 

Aluminum 



Property 

Value 

Comment 

1 

1.0 

Vacuum gap 

2 

0.001 

'vlO' 2 mesh units 

3 

-1. 

Conductor 

4 

13. 


5 

0.97 ) 


6 

0.3 keV ( 

See Reference 2-4 

7 

260 A ^ 


8 

!.3 

See Reference 5 

9 

240 A ) 


10 

1.73 


11 

1.76 


12 

40 keV 


13 

4.0 x 10" 5 A/m 2 

See Reference 6 
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Magnesium 

Property 

1 

2 

3 

4- 

5 

6 

7 

8 
9 

10 

11 

12 

13 


Value 

1.0 

0.001 

- 1 . 


12 . 

0.92 

0.25 keV 

- 1 . 


1.74 
24.3 
1.36 
40. keV 

2.0 x 10~ 5 A/m 2 


Comment 

Vacuum gap 
~2 

'v-lO mesh units 
Conductor 


See Reference 2-4 

Feldman's formulation, 
Reference 1 

Density, gm/cm" 

Atomic weight 


Teflon 

Property 

1 

2 

3 

4 

5 

6 

7 

8 
9 

10 

11 

12 

13 


value 

2.0 

1.27 x 10~ 4 
1.0 x 10 ^ 4 mhc/m 


10 . 


3.0 ) 

0.3 keV i 


- 1 . 


2,00 
16.7 
1.40 
70 keV 

2.0 x 10~ 5 A/m 2 


Comment 


5 mil 


Reference 7 

Feldman's formulation , 
Reference 1 

3 

Density, gm/cm 
Atomic weight 
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Kapton 


Property 

Value 

Comment 

1 

3.5 


2 

1.27 x 10“ 4 

5 mil 

3 

1.0 x 10“ 14 mho/m 


4 

5.0 


5 

2.1 \ 

Reference 7 

6 

0.150 keV J 


7 

8 

-1. 

Feldman's formulation 
Reference 1 

"***■ 


9 

1.42 

Density, gm/cm 2 

10 

9.8 

Atomic weight 

11 

1.40 


12 

70. keV 


13 

2.0 x 10“ 5 A/m 2 


sio 2 



Property 

Value 

Coiranen c 

1 

4.0 


2 

1.27 x 10~ 4 

5 mil 

3 

1.0 x 10 ~ 1 4 mho/m 


4 

10. 


5 

2.4 | 

0.4 keV | 

Reference 2-4 

6 


7 

250 A } 


8 

0.12 [ 

Reference 8 

9 

360 A ' 


10 

1.63 


11 

1.40 


12 

70.0 keV 


13 

2.0 x 10~ 5 A/m 2 
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6.6 CODE DETAILS 

The main routine for the object definition section of 
NASCAP is the subroutine OBJDEF. A flowchart of this routine 
is shown in Figure 6.18. The three main functions are: 

(1) determine the geometrical configuration and surface cell 
list from user input, (2) perform preprocessing of the material 
properties, and (3) determine the surface points, special 
points and their potential coefficients. The subroutine 
writeups which follow are divided into (1) general and mis- 
cellaneous routines, (2) geometrical object routines, (3) 
material routines, and (4) point and potential coefficient 
routines . 

6.6.1 General and Miscellaneous Routines 
CMPRSS 

Purpose: Eliminate superfluous or redundant surface cells. 

Calling Sequenc e: SUBROUTINE CMPRSS ( LTBL , NX , NY , N Z ) 

Called By : INPUT, OBJDEF, QSPHER 


CONDUC 

Purpose : Places conductor indices in surface cell array. 

Calling Sequence : SUBROUTINE CONDUC (NCON) 

Called By : INPUT 


DIAGNO 

Purpose : Prints list of surface cells. 

Callincr Seauence: SUBROUTINE DIAGNO 


Called By : OBJDEF 
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Fiyure 6.10. OBJDEF flow chart 



















DELETE 


Purpose ; Deletes a volume cell and all associated surface 
cells. 

Calling Sequence : SUBROUTINE DELETE (IOFF ,LTBL,NX,NY,NZ , IJK) 

Called By : INPUT, FIL111 


INPUT 

Purpose : Reads keywords from ISAT input file and calls objec 

definition routines. 

Calling Sequence ; SUBROUTINE INPUT (IOFF, LTBL, NX , NY ,NZ , ISAT) 
Called By ; OBJDEF 


LORDEQ 

Purpose : Identical to LORDER (q.v.). 

Calling Sequence ; SUBROUTINE LORDEQ (LIST,W) 
Called Bv: OBJDEF 


NUMLTB 

Purpose : Number volume cells (LTBL) adjacent to ob]ect for 

purpose of calculating low-energy electron flux. 

Calling Sequence : SUBROUTINE NUMLTB (LTBL, NX ,NY ,NZ) 

Called By: OBJDEF 
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OBJDEF 


Purpose : Primary routine for object definition section of 

NASCAP code. (See Figure 6.18.) 

Calling Sequence : SUBROUTINE OBJDEF (NX ,NY,NZ ,NC) 

Called By : NASCAP 


SETMC 

Purpose : Establishes the surface point as the J tJl -multi- 

conductor point. 

Calling Sequence : SUBROUTINE SETMC (I,J) 

Called By : SRFPTS 

6.6.2 Geometrical Object Routines 
CFLICT 

Purpose : Determines if a right-triangular surface cell is 

valid, superfluous or illegal . 

Calling Sequence ; SUBROUTINE CFLICT (N100 ,NORM , IB5 , T.RET) 

Called By ; Several geometrical object routines and CMPRSS. 

CUBE34 

Purpose : Defines a truncated-cube volume cell. 

Calling Sequence : SUBROUTINE CUBE 3 4 (IOFF , LTBL,NX ,NY ,NZ , STUFF ) 

Called By: FIL111 
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FILlll 


Purpose : Defines a wedge with (111) -type surface normal. 

Calling Sequence : SUBROUTINE FILlll (IOFF ,LTBL,NX ,NY , NZ , ISAT ) 

Called By : INPUT 

NIOOBJ 

Purpose : Defines a rectangular parallelepiped with parameters 

passed in array JSTUFF. 

Calling Sequence : SUBROUTINE NIOOBJ (LTBL,NX , NY ,NZ , JSTUFJJ^ 

Called Bv: NIOOCT 


NIOOCT 

Purpose : Defines an octagonal right cylinder as two RECTAN'l 

and four WEDGE's from parameters passed in array JSTUFF. 

Calling Sequence : SUBROUTINE NIOOCT (IOFF , LTBL, NX , NY, NZ , JSTUFF ) 

Called By : OCTGON, QSPHER 


NIOTET 

Purpose : Defines a tetrahedron from parameters passed in 

array JSTUFF. 

Calling Sequence : SUBROUTINE NIOTET ( IOFF, LTBL, NX, NY, NZ , JSTUFF) 

Called B y: QSPHER 
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NIOWGE 


Purpose : Defines a right triangular cylinder (WEDGE) from 

parameters passed in array JSTUFF. 

Calling Sequence : SUBROUTINE NIOWGE (LTBL, NX, NY, NZ , JSTUFF) 

Called By : OCTGON, QSPHER 


OCTGON 

Purpose : Defines an octagonal right cylinder from input data 

as two rectangular parallelepipeds and four wedges. 

Calling Sequence : SUBROUTINE OCTGON (IOFF,LTBL,NX,NY,NZ ,ISAT) 

Called By : INPUT 

QSPHER 

Purpose ; Defines a quasi-sphere from input data. 

Calling Sequence ; SUBROUTINE QSPHER (IGFF , LTBL,NX,NY,NZ , ISAT) 
Called By : INPUT 


RFC TAN 

Purpose : Defines a rectangular parallelepiped from input data. 

Calling Sequence : SUBROUTINE RECTAN (IOFF # LTBL, NX , NY , NZ , ISAT) 

Called By : INPUT 


TETRAH 

Purpose : Defines a tetrahedron from input data. 

Calling Sequence; SUBROUTINE TETRAH (IOP F , LTBL, NX ,NY , NZ , ISAT) 


Called By : INPUT 
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TRNGLS 


Purpose : Performs further processing of right-triangle surface 

cells / including determination of right-angle corner locatirn 
and elimination of redundant cells. 

Calling Sequence : SUBROUTINE TRNGLS (LTBL,NX,NY,NZ) 

Called By ; OBJDEF 


WEDGE 

Purpose : Defines a right triangular cylinder (WEDGE) from 

input data. 

Calling Sequence : SUBROUTINE WEDGE (IOFF ,LTBL, NX , NY ,NZ , ISAT) 

Called By : INPUT 

6.6.3 Materials Routines 
ELSEC1 

Purpose : Identical to ELSEC (q.v.) 

Calling Sequence i SUBROUTINE ELSECl (YIELD , E ,THETA, IOPT) 
Called By : INEIMP 


INEIMP 

Purpose ; Preprocesses material parameters associated with 
electron secondary emission. 

Calling Sequence : SUBROUTINE INEIMP (MPRM) 

Called By; MATPRO 


165 


i 



MAT PRO 


Purpose ; Preprocessing of material property parameters. 
Calling Sequence : SUBROUTINE MATPRO 

Called By ; OBJDEF 


QEQN 

Purpose : Routine used in preprocessing of electron-produced 

secondary material parameters. 

Calling Sequence : FUNCTION QEQN (Z,K,Q) 

Called 3v : INEIMP through ZSYSTM 

QEQN 2 

Purpose : Routine used in preprocessing of electron-produced 

Secondary material parameters. 

Calling Sequence : FUNCTION QEQN 2 (Z,K,Q) 

Called Bv : INEIMP through ZSYSTM 


ZSYSTM 

Purpose : Modified version of IMSL routine to solve a system 

of nonlinear equations. 

Calling Sequence : SUBROUTINE ZSYSTM (F , EPS ,NSIG ,N,X , ITMAX , 

WA , PAR , IER) 

Called Bv: INEIMP 
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6.6.4 Point and Potential Coefficient Routines 
ISIT 

Purpose : Extracts N bits from word I, with bit NO the right- 

most bit. (The rightmost bit of a word is bit O.) Substitute 
for FLD functions for machine independence. 

Calling Sequence : FUNCTION IBIT (I,NO,N) 

Called By: PERMU , SPECEL, others 


PERMU 


Purpose ; Generates permutation vector from LT3L entry. 
(See Figure 6.15.) 

Calling Sequence : SUBROUTINE PERMU (IPMU,ECODE) 

" IPMU (8 ) - Permutation vector for 

cube vertices 

ECODE - element table (LT3L) entry 


Called Bv: SPCWGT 


SPECEL 

Purpose ; Identifies special cells and puts LTBL m its final 
form except for numbering by NUMLTS. (See Figures 6.19-5.24.) 

Calling Sequence ; SUBROUTINE SPECEL' (LTBL, NX ,NY ,NZ ) 

Called Bv: OBJDEF 


SPCPTS 

Purpose : Forms list of special points. 

Calling Sequence : SUBROUTINE SPCPTS (LTBL, NX ,NY,NZ) 

Called Bv: OBJDEF 
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(Format) 

Description 
Standard Orientation 

Potential Function = 

Weight Matrix, ^ = £ ^ij^i^j 


Point Index 

1 

2 

3 

4 

5 

6 

7 

8 


Cube Corner 

0 0 0 
10 0 
0 10 
110 
0 0 1 
10 1 
0 11 
111 


Figure 6.19. Format and definition of point indices for 
Figures 6.20-6.24. 
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Standard Cell 0 


Empty trilinear cube 
Orientation: Arbitrary 


Potential Function: 


1 
1 

2 

3 

4 

5 

6 

7 

8 


(1-x) (1-y) ( 1 - 2 ) 
(1-z) (1-y) x 
(l-x)y (1-z) 
(1-z) yx 
2 (r-y) ( 1— x) 
x(l-y ) (2 
zy (1-x) 
xyz 



7 




1/3 


0 

0 

-1/12 

0 

1/3 

-1/12 

0 

-1/12 

1/3 

0 

-1/12 

1/3 

-1/12 

1/3 


-1/12 

0 

-1/12 

-1/12 

0 

1/3 

-1/12 

-1/12 

0 

-1/12 

0 

-1/12 1/3 

-1/12 

-1/12 

-1/12 

0 

-1/12 

0 0 


Ficrure 6.20. The empty cube volume cell. 
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Standard Cell 1 

Hair-Empty Wedge 

1 < x + y < 2 
0 < 2 < 1 

Orientation: Right angle along 

line 7-8 


Potential function: 

1 N 1 

1 0 

2 (1-y) (1-z) 

3 Cl-x) (1-z) 

4 (x+y-1) (1-z) 

5 0 

6 (1-y) z 

7 (l-x)z 

8 (x+y-1) z 



0 


0 

1/4 


0 

1/24 

1/4 

0 

-1/8 

-1/8 

0 

0 

0 

0 

0 

-1/24 

0 

-1/24 

0 

0_ 

-1/8 

-1/8 


5/12 



0 

0 


-1/8 

0 

1/4 

-1/8 

0 

1/24 

-1/12 

0 

-1/8 


1/4 

- 1/8 


5/12 


Figure 6.21. The half-empty wedge volume cell. 
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Standard Cell 2 

Cube with diagonal line on one face 
Orientation: Line from 2 to 3 



Potential Function: 
i N 1 

1 (1-x-y) ( 1 - 2 ) 0 (1-x-y) 

2 fx0 (1-x-y) + <l-y) 0 (x+y-1) 3 (1-z) 

3 [y0 (1-x-y) + (1-x) 0 (x-*-y-l) ] (1-z) 

4 (x+y-1) (1-z) 0 (x+y-1) 

5 (1-x) (1-y) z 

6 x(l-y)z 

7 (l-x)yz 

8 xyz 



5/12 

- 1/8 1/2 
- 1/8 1/12 
0 - 1/8 
7/360 -37/360 

-11/180 -1/45 
-11/180 -19/180 
-23/360 -37/360 


1/2 

- 1/8 5/12 
-37/360 -23/350 1/3 
-19/180 -31/180 0 
-1/45 -11/180 0 

-37/360 7/360 -1/12 


1/3 

- 1/12 

0 


1/3 

0 


1/3 


Figure 6.22. The empty special cell. 


171 


i 


I 




Standard Cell 3 

Tetrahedron 

2<x + y+ z<3 

Orientation: Empty corner at 

point 8 

Potential Function: 

1 N 1 

1 0 

2 0 

3 0 

4 1-z 

5 0 

6 1-y 

7 1-x 

8 x+y+z-2 


0 

0 0 

0 0 1/6 

0 0 0 0 

0 0 0 0 1/S 

0 0 0 0 0 1/6 

0 0 - 1/6 0 - 1/6 - 1/6 1/2 

Figure 6.23. The tetrahedron volume cell. 
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Standard Cell 4 


Truncated Cube 

Orientation: 000 corner (point 1) 

missing 



Potential Function: 


1 
1 

2 

3 

4 

5 

6 

7 

8 


N L 

0 

exercise 

for 

reader 



0 


0 

5/12 


0 

1/72 

5/12 

0 

-11/120 

-37/360 

0 

1/72 

1/72 

0 

-37/360 

-1/9 

0 

-1/9 

-11/120 

0 

-5/36 

-5/36 


13/36 


-1/9 

5/12 

-7/180 

-11/120 13/36 

-7/180 

-37/360 -7/180 

1/4 5 

-5/36 1/45 


13/36 

1/45 


Figure 6.24. The empty truncated cube volume cell. 
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SPCWGT 


Purpose : Calculates potential coefficients for surface and 

special points associated with volume cells. (See Figures 
6 . 19 - 6 . 24 .) 

Calling Sequence : SUBROUTINE SPCWGT (LTBL,NX,NY,NZ ) 

Called Bv: OBJDEF 


SRFPTS 

Purpose : Identifies surface points and computes potential 

coefficients associated with dielectric layers and vacuum gaps. 
Identifies multiconductor points. 

Calling Sequence : SUBROUTINE SRFPTS (LTBL,NX ,NY ,NZ ) 

Called Bv: OBJDEF 



7. CODE DESCRIPTION - POTENTIAL SECTION 


7.1 CONJUGATE GRADIENT POTENTIAL SOLUTION 

Potentials at grid points and on each of the conduc- 
tors are solved using a conjugate gradient (CG) iteration 
scheme. *■ J ' Either of two slightly different solution 

algorithms nay be specified by the user. The first, referred 
to as algorithm 0, is slightly faster than the second and re- 
quires one less external storage file, bur is less convenient-, 
to use because there is no straightforward way to determine 
quality of convergence. A generally decreasing trend in the 
residual may often be mashed by large short-tem oscillations 
over several orders of magnitude. The second algorithm, how- 
ever, referred to as algorithm 1, constructs a monotonically 
decreasing residual quantity which allows the user to more 
accurately judge the quality of the' potential convergence- at 
any point in the potential solution. In both cases, the 
residuals and ether relevant quantities are printed at each 
potential iteration and also printer-plotted over all itera- 
tions at the end of the potential iteration loop for ease of 
examination. 

Although each actual application will demand specific 
judgments, the following guidelines will help the user m 
deciding which algorithm to use: (1) Algorithm I is recom- 

mended for most applications. Where a new or unfamiliar case 
is being investigated, the mor.ctcnic residual is a great ccr.ve 
nience in judging accuracy of the potentials. The extra com- 
puting time involved in this algorithm is small. '2) In cer- 
tain cases where many runs have been made on a given case or 
on similar cases, and where user judgment as to how many 
iterations will produce usable potentials is fairly solid, 
algorithm 0 would be a proper choice if costly large-scale 
production runs are to be made. In this case, algorithm l’s 
extra overhead in maintaining the additional external storage 
file and communicating with it may add a small (percentage- 



wise) but significant amount of dollars to the total cost. 

In addition, although algorithm 1 yields a monotonic residual 
and algorithm 0 may yield an oscillating but generally de- 
creasing residual, the actual quality of convergence of 
algorithm 0 at any given potential iteration will quite 
likely be the same as that derived from algorithm 1. The 
difference is simply that the user obtains a clearer picture 
of the quality of convergence with algorithm 1. When any 
uncertainty about the convergence of the potentials exists or 
is expected, algorithm 1 should definitely be used. 

Requirements of the conjugate gradient solution tech- 
nique are a system of linear equations whose coefficient array 
is both symmetric and positive definite. Our 3-D Poisson 
equation system is arranged to meet these requirements. Given 
a system of N unknowns, the conjugate gradient equations will 
converge to an exact solution in N iterations and will 
apparently tend to converge in most cases to a practical solu- 
tion after a few times iterations. 

The solution equations for algorithm 0 are as fol- 
lows : 


Given 


Ap = ROUS 


(7.1) 


where A is a coefficient array, p is the potential (solution) 
vector, and RCUS is a boundary condition vector; we need to 
maintain storage for vectors (Ap) , p, ROUS, two more vectors, 
r and u, and two scalars, a and 2. Solution proceeds as follows 


Initialization: 


u° = ROUS - ( Ap ) ° 


o o 
r = u 


(7.2) 

(7.3) 


n = 0 


(7.4) 
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(7.5) 


Iteration Loop: (Au) n = (A) (u n ) 


n _ 

r *r 

(7.6) 

Ct — 

u°‘ (Au) n 

n+1 

p 

= p n + ot n u n 

(7.7) 

r n+1 

n n , > v i*i 

= r - a (Au) 

(7.8) 

B n = 

-r* (Au) n /u n • (Au) n 

(7.9) 

n+1 

n+1 , „n n 

(7.10) 

u 

= r + 6 u 

n = n 

. + 1 

(7.11) 


Repeat iteration loop 

Figure 7.1 illustrates in block diagram form a gross 
overview of the routines called from subroutine POTENT 
(potential solver routine) which accomplish the conjugate 
gradient solution for algorithm 0. 

The coefficient array A is never explicitly formed,; only 

o n 

elements of the coefficient product arrays (Ap) or (Au) are 
ever formed and stored. Even so, we must form, store, and up- 
date fo~ each iteration five vectors with length equal to the 
number of unknowns in the system. Assuming a moderately fine 
cell spacing of 17 x 17 x 33 points per end, and assuming a 
nesting of three or four of these grids, it is clear that we 
are dealing with significantly large storage reauirements 
(3 x 17 x 17 x 33 = 47,685 vector elements for the innermost 
grid alone , not counting up to seven additional unknowns for 
the conductor potentials) . Obviously, use of the conjugate 
gradient scheme involves a nontrivial data management system 
in order to treat a problem of this scope efficiently and 
still execute in t65K on a UN I VAC 1103 (cr 1110) computer sys- 
tem. 
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To minimize core storage requirements we require tnat 
as few as possible of the five necessary arrays be in core a- 
any given time, and that as few grids worth of data for any of 
those five quantities be in core at any given time. (One 
"grid's worth" of data is NX*NY*N2 numbers which are the 
values of one of the five quantities at all points of any 
one of the several nested NX x NY x NZ grids.) 

Five external files (drum, disc or tape) are utilized 
for storage of the five large vector quantities; one additional 
file acts as a "rotating spare" for use in the efficient up- 
dating of the p, r and u arrays each iteration. These files 
are designated by variables IP, IR, IU, IAUN, IROUS and 
ISPARE which contain the external file unit numbers. In all 
cases, the quantities are stored on their respective files 
beginning with inner grid data and working outward through the 
nested grids to the outermost grid data at the end of the file. 
All data transfer between these files and the code buffer ar- 
rays is done by a subroutine called I10VDAT. MOVDAT transfers 
data both into core from external files and vice versa, and 
in all cases uses efficient block transfer routines to move 
an entire grid's worth of data at one time. 

Figure 7.2 illustrates the sectioning of the nested 
grids for two adjacent layers. The limiting case for simulta- 
neous core storage of grid data occurs m subroutine COPROD 
which, as shown in Figure 7.2, calculates" either the (Apl° ar- 
ray (initialization) or the (Au) n array (within iteration loop) . 
For purposes of explanation, we will speak of the resultant 
coefficient product array here simply as (Au) *. A self- 
consistent starting point at an arbitrary grid m the middle 
of the next would assume I-A and all points within it have 
already been dealt with. Calculation of Area I-B of (Au) r ‘ 
requires only two grids worth of simultaneous storage: one 

array containing . u data for that grid (remember that the 
A matrix is never formed explicitly — its elements are 
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Figure 


A: Inner edge of a grid 

B: Interior space of a grid 

C: Outer edge of a grid 

I: Inner grid of the two 

0: Outer grid of the two 



Key: 

— Grid boundaries 

Area I-A Inner Edge of Inner Grid 


\\\\\ 


///// 


Area I-B Interior Space of Inner Grid 

Area I-C Outer Edge of Inner Grid 
Area 0-A Inner Edge of Outer Grid 

Area 0-B Interior Space of Outer Grid 

Area 0-C Outer Edge of Outer Grid 


7.2. Simplified two-dimensional sectioning of two 
adjacent grids. 
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essentially calculated one at a time as needed) , and a second 
array into which are placed the resultant elements of (Au) n 
for that grid. However, the calculation of elements of (Au) n 
for either Area I-C or Area 0-A requires u data for both grids 
to be present simultaneously. We place these data in buffer 
arrays BUF1 and 3UF2. Fortunately, we need only maintain one 
grid's worth of resultant storage for (Au) n at any one time, 
since once Areas I-A, I-B and I-C of (Au) n have been calcu- 
lated, then BUF3, the array containing the resultant (Au) n , 
can be written out onto its external storage file by MGVDAT. 
Then, with the same two grids' worth of. u data in BUF1 and 
BUF2, the Area 0-A elements of (Au) n may be calculated and 
placed in BUF3 , overwriting the last grid’s data which has 
just been written out. The whole scheme, given several 
grids in the problem, is essentially a sort of piggy-back- 
leapfrog data manipulation which might be diagrammed as m 
Figure 7.3. 

Thus we see that the mcst data we need ever have m 
core at any one time is three grids' worth; two of u (or p 
in the initialization) , and one of Au n (or Ap°) . The synop- 
sis of subroutine COPRQD presented m Figure 7.4 illustrates 
the actual mechanism of implementing what is diagrammed m 
Figure 7.3 including all data transfers between cere and 
external storage files. 

All other routines used by the conjugate gradient 
potential solver require at most two grid's worth of data at 
any one time. A synopsis of PUPDAT, a representative routine 
of this group, is presented in Figure 7.5. 

Subroutines RUPDAT, UUPDAT and URSETO follow essentially 
the same pattern as shown here. The only concept of note m 
the three update routines is the use of the "rotating spare" 
file. Eor each of PUPDAT, RUPDAT, and UUPDAT, a quantity is 
being used to update itself: e.g., p n+ ^ = p“ + a n u n . Re- 

call that despite the simplicity of this equation, each cf 
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Grid 

Index 


U-Data 

Blocks 


Grid 

Pair 


1 


2 


3 


4 


5 



Calculate (Au) n for 1-B 
(Area A is not meaningful for 
the innermost grid.) 


ilculate 
rite out 


Calculate 
Write out 
Calculate 


Calculate 
Write out 
Calculate 

Calculate 
Write out 
Calculate 

Write out 


(Au) n 

fer 1-C 

(Au) n 

to unit IAUN 

<Au) n 

for 2-A, 2-B 

(Au) n 

for 2-C 

(Au) n 

to unit IAUN 

(Au) n 

for 3-A, 3-B 

(Au) n 

for 3-C 

(Au ) n 

to unit IAUN 

(Au)n 

for 4 -A, 4-3 

(Au) n 

for 4-C 

(Au ) n 

to unit IAUN 

(Au) n 

for 5-A, 5-B 

(Au ) n 

to unit IAUN 


(Area 5-C, the outer edge of the outer lost grid, is a specia 
case treated elsewhere by the cede and need not be filled in 
here. ) 


Figure 7.3. 


Data manipulation in the COPROD routine m a 
five grid problem. 
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Subroutine COPROD(IA) 


4 


- $- 

« i 

, 3 

t 


... <f- 


V 


■f 


IA is an I/O unit number containing the file identi- 
fier for either the p array (initialization) or the u array 
(iteration loop) . The description here assures we are 
dealing with the u array for this call. 


Call MOVDAT: Innermost grid u data from file 


IU 


BUF1 


Call ISPACE: 


, n 


{ Au) “ for Area 1-B - BUF3 
(Area 1-A has no meaning.) 


Loop over nested grids: DO 1G=2,NG 


—Ml = 1 

N2 = NX *NY*NZ+1 

If loop index 1G is odd, switch values of XI and N2 


Call MOVDAT: 
Call IRC ALL: 

Cali MOVDAT: 

Call ORCALL: 


Next grids u data -* BUF1(N2) 

Calculate Area. I-C of (Au) n using u 
data for two adjacent grids in EUFi(Nl) 
and BUF1 (N2) ; results go m 3UF3 

Write completed inner grid (Au) n dat3 
(Areas I-A, I-B and I-C) from 
BUF3 -* File IAUN 

Calculate Area O-A of (Au) n using u 
data for two adjacent grids; results 
go in BUF3 , overwriting previous con- 
tents 


Repeat loop for next grid pair 


Call MOVDAJ 


Write outermost grid (Au 
BUF3 - File IAUN 


; n 


cata zrcn 


Rewind IA (IP or 
Rewind IAUN 
Return 


:u) 


Figure 7.4. Synopsis af subrcutme COPROD. 
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PUPDAT does the equation p n+ ^ = p n + a n u 11 


-Loop over nested grids: Loop index IG=1, NG 

Call MOVDAT: Next or ids p,IP -*■ BUF1 
Call MOVDAT: Next grids u,IU - BUF2 

Loop over each point in the grid 

i— *-BUFl (I,J,K) = BUF1 {I, J, K) + ALPHA*BUF2 (I, J,K) 


All points' 
. done? ^ 


MOVDAT: BUF1 (p ] -*■ ISPARE 


No Yes 

-^ Done grids?^> 


IP ISAVE 
ISPARE -*■ IP 
ISAVE -* ISPARE 

Rewind IP, IU, ISPARE 

Return 



Synopsis of subroutine PUPDAT (representative 
of RUPDAT and UUPDAT as well) . 
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p and u is actually a vector of perhaps 38,148 elements (for 
a four grid problem) , treated in blocks of 9537 elements at a 
time (one grid's worth for a 17 x 17 x 33 grid). As is il- 
lustrated in Figure 7.5, one grid's worth of p elements and 
one of r elements are brought into arrays BUF1 and BUF2. The 
updated p n+ ^ elements for that grid are calculated and then 
must be written out onto some external storage file. (Writing 
them into the unused BUF3 array would not eliminate the problem, 
which is one of preserving proper ordering and spacing of data 
on the external storage files.) The data obviously cannot be 
written back into the p file from which p n was just read, or 
the next grid's worth of data would be overwritten. One could 
devise a cumbersome rewind and space-forward-to-f lag technique, 
but the coding work and, more significantly, data read-write 
time involved in such an approach would not be at all justi- 

n+i 

fied given the following alternative: The updated p is 

written out onto whatever unit is designated as unit ISPARE 
at that time. The data is written out in the same natural 
order in which it was read in and calculated, one grid at a 
time, beginning witn the inner grid data. When all computa- 
tions are finished, we have old file IP at the end of its p n 
data, file IU at the end of its u n data, and file ISPARE at. 
the end of its updated p n+1 data. To recreate the natural 
situation expected by the rest of the code — new p data on 
IP, u data on IU, and ISPARE a free file with no longer 
needed data on it, we need only switch the values held by 
variables IP and ISPARE. Since the values of IP, ISPARE, IR, 

IU, IAUN and IB are known to the rest of the code, and since 
the old p n values on the new file ISPARE are indeed no longer 
needed for any calculations, we can do this with impunity. 

All that remains is to rewind the three files IP, IU and 
ISPARE so that at the start of any subsequent routine, we 
are situated at the opening (innermost grid) data for each 
file, as is required by all routines in the code. 
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We see, then, that the four I/O unit variables IP, IR, 
10 and ISPARE do, in fact, "rotate around" the set of identi- 
fiers for the corresponding external files through the course 
of several iterations, and no ambiguity ever exists as to 
which file contains which data set. In this way, the trivial 
sacrifice of maintaining one more external file during the 
course of a run allows simple and efficient coding for and 
execution of the three array updates for p, r and u. 

The execution of algorithm 1 is quite similar, but in- 
volves one more vector which we call y, stored on unit IY, 
and an additional scalar, c. 

The equations for algorithm 1 are as follows: 


Initialization: 

u° = ROUS - (Ap)° 

(7.12) 


o o 

r = u 

(7.13) 


o .o 

y = <Ap) 

(7.14) 


H 

n 

o 

o 

(7.15) 

Iteration Loop: 

(Au) r - = (A) (U n ) 

(7.16) 


n r n *r n 

a = — 

n , > v n 
u ° (Au) 

(7.17) 


n+l n . n n 

P = p + u p 

(7.18) 


n+l n n . n 

r = r - a (Au) 

(7.19) 


e n = -r ■ (Au.t n /u" ■ (Au) n 

(7.20) 


c n+1 = 1 + 0 n c n 

(7.21) 
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Calculate Monotonic , 

Residual: x*r = u*u n /c (7.22) 

y n+1 = ROUS - r n+1 + S n y n (7.23) 
u n+1 = c n+1 ROUS - y n+1 (7.24) 

n = n + 1 

Repeat Iteration-Loop 

Figure 7.6 illustrates in block diagram form the execu- 
tion of algorithm 1. Aside from the calculation of the mono- 
tonic residual (which has no effect on the solution vector) 
the only differences in execution between the two algorithms 
are the initialization and updating of the y vector and c 
scalar, and a change in the updating of the u vector. Sub- 
routines YINIT and YUPDAT perform the y vector operations. 
Subroutine UUPDAT updates the u vector for algorithm 1, while 
subroutine UUPDAU updates the u vector for algorithm 0. In 
addition, file IY now rotates unit numbers along with IP, IR, 

IU and I SPARE. 


Following is a complete documented list of all sub- 
routines involved in the potential solution. Subroutine 
arguments not defined under the routine in which they appear 
are general quantities which are defined in the variable 
glossary. 


187 
































I 

\ 


\ 


V. 


I 




7.2 SUBROUTINE DESCRIPTIONS 
COPROD 

Calling Sequence : SUBROUTINE COPROD (IA, IAUN, BUF1, BUF3 , 

NG, NX, NY, NZ, P3, W4 , UDOTAU, AUNCON, 
UPCOND, IOBJ, H, LIST, WPCOND, PMCND, 

SWPSSC, CIJ, CIJSUM, MBIAS, MFIX) 

Called By : POTENT (in two places) 

Local Variables : IA: COPROD takes the product of the coef- 

ficient array A with whatever vector 
is on external storage unit IA. In 
the initial call to COPROD, IA = IP, 
the unit storing the (at that time) 
initial guess for the potentials. In 
the iteration loop call to COPROD, 

IA = IU, the unit storing the u n vector. 

Purpose : COPROD computes the product of the linear equation 

coefficient array A with either of p° or u n , depending on 
whether the call to COPROD is the initial call in POTENT or 
the call within the potential ite: ation loop. COPROD also 
calculates u n °(Au) n adding terms into the dot product one at 
a time, as it calculates the individual terms of (Au) n . 

Internal Organization : COPROD moves the inner grid's worth 

of data from IA into BUF1, and calls ISPACE and OBJSPC to com- 
pute the coefficient product values at grid points within the 
1— B (see Figure 7.2) grid area (placed in BUF3 ) . 

A loop over the number of nested grids follows. This 
loop reads in the next grid's IA datJ , calls IRCALL to calcu- 
late product values at grid points within area 1-C (outer 
edge of inner grid) , writes the completed inner grid (Au~) n 
data from BUF3 onto file IAUN, and finally, calls ORCALL to 
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calculate product values for points within grid areas 2-A 
and 2-B (inner edge of outer grid and interior space of 
outer grid) after BUF3 has been re-initialized to zero) . 

The loop repeats itself for all grids, at each point having 
two adjacent grid's IA data plus one grid's IAUN data in 
core at one time, and making its way one grid at a time 
outward towards the outermost grid pair. This procedure 
is diagrammed in Figure 7.3. 

Finally, both the IA and IAUN files are rewound. 
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ICORN 

Calling Sequence : SUBROUTINE ICORN (PI, P2, BUF3 , NX, NY, 

HZ, P2D, W4, UDOTAU) 

Called By ; IRCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space BUF2 

P2D: Storage space P3 (scratch storage) 

Purpose : Using the trilinear incerpolation scheme, ICORN 

computes the coefficient product values for the grid points 
of the outer edge- of the inner grid of the current pair (grid 
area I-C) which are the eight outer corner points of the inner 
grid. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3. 
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INEDGX 

Calling Sequence : SUBROUTINE INEDGX (PI, P2, BUF3 , NX, NY, 

NZ, P2D, W4 , UDOTAU) 

Called By ; IRCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space BUF2 

P2D: Storage space P3 (scratch storage) 

Purpose : Using the trilinear interpolation scheme, INEDGX 

computes coefficient product values for the grid points of 
the outer edge of the inner grid of the current pair (grid 
area I-C) which lie on edges (where two planer meet) parallel 
to the x-axis. Corner points on these edges are not treated 
here. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u 11 (or p°) . Coefficient product 
values are placed in BUF3 . 
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Calling Sequence : SUBROUTINE INEDGY (PI , P2, BUF3 , NX, NY, 

NZ, P2D, W4 , UDOTAU) 

Called By : IRCALL 

Local Variables ; Pis Storage space BUFl 

P2: Storage space BUF2 

P2D: Storage space P3 (scratch storage) 

Purpose : Using the trilinsar interpolation scheme, INEDGY 

computes coefficient product values for the grid points of 
the outer edge of the inner grid of the current pair (grid 
area I-C) which lie on edges (where two planes meet) parallel 
to the y-axis. Corner points on these edges are not treated 
here. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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INEDGZ 

Calling Sequence : SUBROUTINE INEDGZ (PI, P2, BUF3 , NX, NY, 

NZ, P2D, W4 , UDOTAU) 

Called By ; IRCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space BUF2 

P2D: Storage space F3 (scratch storage) 

Purpose : Using the trilinear interpolation scheme, INEDGZ 

computes coefficient product values for the grid points of 
the outer edge of the inner grid of the current pair (grid 
area I-C) which lie on edges (where two planes meet) parallel 
to the z-axis. Corner points on these edges are not treated 
here. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3. 
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IPLANX 

Calling Sequence : SUBROUTINE IPLANX (PI, P2, BUF3 , NX, NY, 

NZ, P2D, W4 , UDOTAU) 

Called By : IRCALL 

Local Variable s; PI: Storage space EUF1 

P2: Storage space DUF2 

P2D: Storage space P'i (scratch storage) 

Purpose : Using the trilinear interpolation scheme, IPLANX 

computes coefficient product values for two plane sections 
of the outer edge of the inner grid of the inner pairs (grid 
area I-C) : (1) the plane of grid points on the +x end of 

the inner grid, except for the edge points and corner points, 
and (2) the plane of grid points on the -x end of the inner 
grid, except for the edge points and corner points. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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IPLANY 

Calling Sequence ; SUBROUTINE IPLAIJY (PI, P2, BUF3 , NX, NY, 

NZ, P2D, W4 , UDOTAU) 

Called By ; IRCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space 3UF2 

P2D: Storage space P3 (scratch storage) 

Purpose : Using the trilinear interpolation scheme, IPLANY 

computes coefficient product values for two plane sections 
of the outer edge of the inner grid of the current pairs 
(grid area I-C) : (1) the plane of grid points on the +y end 

of the inner grid, except for the edge and corner points, and 
(2) the plane of grid points on the -y end of the inner grid 
except for the edge points and corner points. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°). Coefficient product 
values are placed i.. BUF3 . 
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IPLANZ 

Calling Sequence : SUBROUTINE IPLANZ (PI, P2, BUF3, NX, NY, 

NZ, P2D, W4 , UDOTAU) 

Called By : IRCALL 

Local Variab3es : PI: Storage space BUF1 

P2: Storage space BUF2 

P2D: Storage space P3 (scratch storage) 

Purpose : Using the trilinear interpolation scheme, IPLANZ com- 

putes coefficient product values for two plane sections of the 
outer edge of the inner grid of the current pair (grid area 
I-C) : (1) the plane of grid points on the +z end of the 
inner grid, except for the edge and corner points, and (2) 
the plane of grid points in the -z end of the inner grid, 
except for edge and corner points. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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IRCALL 

Calling Sequence ; SUBROUTINE IRCALL (BUF1, BUF 2 , BUF3 , NX, 

NY, N2 , P3 , W4, UDOTAU) 

Called By : COPROD 

Purpose : IRCALL calls the series of I prefix routines which 

compute coefficient product values for grid area I-C in Fig- 
ure 7.2 (outer edge of the inner grid of the current pair). 

The routines called by IRCALL are (in sequence) : 

IPLANX 
I PLAN Y 
I PLAN Z 
INEDGX 
INEDGY 
INEDGZ 
ICORN 
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ISPACE 


Calling Sequence : SUBROUTINE ISPACE (P, BUFJ, NX, NY, N2, 

UDOTAU) 

Called By : COPROD 

Purpose ; ISPACE calculates coefficient product values at grid 
points in the innermost grid of the nested grid system, ex- 
cepting those surface points and special points on or around 
the satellite which we defined in the OBJDEF section and which 
are treated in OBJSPC, and excepting the points on the outer 
edge of the inner grid which are treated by the routines 
called by IRCALL. 



OBJSPC 


Calling Sequence : SUBROUTINE OBJSPC (BUP1, BUF3 , NX, NY, NZ, 

UDOTAU, AUNCON, UPCOND, IOBJ, W, LIST, 
WPCOND, PMCND, SWPSSC, CIJ, CIJSUM, MBIAS, 
MFIX) 


Called By : COPROD 

Purpose : OBJSPC calculates coefficient product values for 

those grid points of the inner grid on or around the satel- 
lite which are defined by the OBJDEF section as being either 
surface points or special points. 

OBJSPC also calculates coefficient product values for 
the unbiased or unfixed conductor potentials. 
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OCORN 

Calling Sequence : SUBROUTINE OCORN (PI, P2, BUF3 , NX, NY, 

NZ, W4, UDOTAU) 

Called By : ORCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space BUF2 

Purpose : Using the trilinear interpolation scheme, OCORN 

computes coefficient product values for grid points at the 
corners of the inner edge of the outer grid of the current 
pair (grid area O-A) . 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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Calling Sequence : SUBROUTINE OPLANX (PI, P2, BUF3 , NX, NY, 

NZ, W4 , UDOTAU) 

Called By ; ORCALL 

Local Variables ; PI: Storage space BUF1 

P2: Storage space BUF2 

Purpose : Using the trilinear interpolation scheme, OPLANX 

computes coefficient product values for two sections of the 
inner edge of the outer grid of the current pair (grid area 
0-A) : (1) the plane of grid points on the +x side adjacent 
to the inner grid of the pair, and (2) the plane of grid 
points in the -x side adjacent to the inner grid of the pair. 
In both cases, corner points and the edge points where two 
planes meet are excepted. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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OPLANY 

Calling Sequence ; SUBROUTINE OPLANY (PI, P2, BUF3 , NX, NY, 

NZ, W4, UDOTAU) 

Called By : ORCALL 

Local Variables : Pis Storage space BUF1 

P2: Storage space BUF2 

Purpose : Using the trilinear interpolation scheme, OPLANY 

computes coefficient product values for two sections of the 
inner edge of the outer grid of the current pair (grid area 
0— A) : (1) the plane of grid points on the +y side adjacent 
to the inner grid of the pair, and (2) the plane of grid 
points on the -y side adjacent to the inner grid of the pair. 
In both cases, corner points and edge points where two planes 
meet are excepted. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values aie placed in BUF3 . 
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OPLANZ 


Calling Sequence : SUBROUTINE OPLANZ (PI, P2, BUF3, NX, NY, 

NZ, W4 , UDOTAU) 

Called By : ORCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space BUF2 

Purpose : Using the trilinear interpolation scheme, OPLANZ 

computes coefficient product values for two sections of the 
inner edge of the outer grid of the current pair (grid area 
0-A) : (1) the plane of grid points on the +z side adjacent 
to the inner grid of the pair, and (2) the plane of grid 
points in che -z side adjacent to the inner grid of the 
pair. In both cacos, corner points and edge points where 
two planes meet are excepted. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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ORCALL 

Calling Sequence : SUBROUTINE OB CALL (BUF1, BUF2 , BUF3 , NX, 

NY, NZ . PJ, W4, UDOTAU) 

Called By ; COPROD 

Purpose : ORCALL first calls the series of O-prefix routines 

which compute coefficient product values for grid area O-A in 
Figure 7.2 (inner edge of outer grid of the current pair), 
and, second, calls subroutine SPACE to compute coefficient 
product values for grid area O-B (interior space of outer 
grid of the current pair) . 

The routines called by ORCALL are (in sequence) : 

OPLANX 

OPLANY 

OPLANZ 

OCORN 

OUTEDG 

SPACE 
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OUTEDG 

Calling Sequence : SUBROUTINE OUTEDG (PI, P2, BUF3 , NX, NY, 

NZ, P2D, W4, UDOTAU) 

Called By : ORCALL 

Local Variables : PI: Storage space BUF1 

P2: Storage space BUF2 

Purpose : Using the tril inear interpolation scheme, OUTEDG 

computes coefficient product values for grid points on the 
inner edge of the outer grid of the current pair (grid area 
0— A) which are on grid edges where two planes meet. Corner 
points of these edges are not treated here. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3. 
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POTENT 


Calling Sequence ; SUBROUTINE POTENT (BUF1, BUF2 , BUF3 , ICY C, 

RITER, UITER/ PCITER, NITERP, ROUS, NX, NY, 
NZ , NG, QCOND, NC, W, LIST, WPCOND, PMCND, 
PCOND, I WARN , IALG) 

Called By : TRILIN 

Purpose : POTENT uses the iterative conjugate gradient tech- 

nique described above to solve a set of linear equations whose 
unknowns are the potentials at all grid points and on each of 
the conductors, and the charges on each biased or fixed- 
potential conductor. 

Internal Organization ; POTENT places fixed-potential conductor 
potentials (set by input in array VFIX) into PCOND, and then 
forms some precalculable quantities (for use in OBJSPC) based 
on the WPCOND and PMCND weights and on the conductor relative 
capacitances. The first call to COPROD calculates A times p° 
where p° is the initial guess for the potentials for this run 
(determined by one of the initial guess options specified by 
user input, or, in the case of a restart, by the final poten- 
tials left over from the previous run) . 

Next, if IALG = 1, YIN IT is called to initialize the 
IY file containing the y vector, and the ROUS vector is 
brought into core from file IROUS. PTROUS is called to cor- 
rect the ROUS array at points in contact with one or more con- 
ductors, and the corrected ROUS is written back out onto file 
IROUS. Next, the ROUSCN array containing ROUS values for the 
conductor points is filled in, and URSETO is called to 
initialize che u° and r° conjugate gradient vectors, and 
to calculate r°*r°. 

The potential iteration loop follows, as described in 
Figure 7.6 above. COPROD is called to calculate A times u n 
for the nth iteration and to compute the dot product u n * (Ac) . 
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The scalar ALPHA is computed from r n *r n and u n *(Au) n , and 
PUPDAT is called to update the solution vector p. Next, 
the r vector is updated by a call to RUPDAT which in addi- 
tion calculates r n+ ^»(Au) n for use in computing BETA. 

If IALG = 1, UDOTUC is called to generate u n *u n for 
use in calculation of the monotonic residual, and YUPDAT is 
called to update the y vector, UUPDAT is then called to per- 
form the update operations on n which are required by algo- 
rithm 1. 

If IALG = 0, UDOTUC, YUPDAT, and UUPDAT are skipped 
and UUPDAU is called instead. UUPDAU performs the u vector 
updating required by algorithm 0. 

All convergence printer plot data for this iteration 
is stored, and the loop repeats itself until HAXITR is 
reached . 

Following completion of the iteration loop (at which 
point a new potential for conductor 1 has been converged 
upon) the potential biases for all biased conductors are 
added to the conductor one potential and placed in their 
respective locations in PCOND. 

Subroutine QCONCP is then called to calculate the 
new amounts of charge on any fixed potential or biased con- 
ductors as determined by the new set of potentials. 

Finally, PTROUS is called again to reset the ROUS 
values to what they were prior to the corrections for points 
on conductors, so that ROUS may be used in its original form 
by subroutine CHARGE and its sequence of routines. 
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PTROUS 


Calling Sequence : SUBROUTINE PTROUS (ROUS, LIST, WPCOND , 

PMCND, VBIAS, IOBJ , MBIAS, NX, NY, IIZ , 
IRESET) 

Called By ; POTENT 

Local Variables ; IRESET: IRESET = 0 implies correct the 

ROUS array for biased conductor 
points. IRESET = 1 implies reset 
the ROUS array to what it was be- 
fore PTROUS corrected it. 

Purpose : The potential equations for all conductors biased 

relative to conductor 1 are added togetj?^5> with the equation 
for conductor 1, terms are collected, and a new single equa- 
tion in one unknown (the potential of conductor 1) replaces 
the original set. 

When this is done the equations for potentials at 
points in contact with biased conductors must also be ad- 
justed. The WPCOND* 9 ^ term in each of these equations be- 
comes WPCONO* (9 + V ix ) = V’JPC0ND*9 1 + WPCOND*V i]L . Since the 
second term is constant, it properly belongs on the right hand 

« 

side of the equation, whose current value is held in ROUS for 
that point. When called with IRESET = 0, PTROUS makes the ap- 
propriate ROUS correction for all points in contact with one 
or more conductors. When called with IRESET = 1, PTROUS re- 
sets ROUS back to what it was before the first -call to PTROUS, 
so that ROUS may be utilized in its original form by sub- 
routine CHARGE, etc. 

It should be noted that PTROUS corrects only the ROUS 
array passed to it; the writing of the corrected values to 
file IROUS is done by POTENT. 
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PUPDAT 


Calling Sequence : SUBROUTINE PUPDAT (ALPHA, IP, IU, ISPARE, 

BUF1, BUF2 , NG, NX, NY, N2, NC, PCOND, 
UCOND, MB I AS, MFIX) 

Called By : POTENT 

Pir pose : PUPDAT updates the linear equation solution vector 

P on file IP (containing system potentials) according to the 
equation 


n+1 _ n 
P = P 


+ a*u 


n 


The updated p n+ *‘ data is written onto file ISPARE, one 
grid at a time, and when all data has been updated, the unit 
numbers IP and ISPARE are switched. 

PUPDAT also updates the conductor potentials in array 

PCOND . 
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QCONCP 

Calling Sequence ; SUBROUTINE QCONCP (QCOND, PCOND, CIJ, 

CIJSUM, WPCOND, PilCND, BUF1, IP, NX, NY, 
NZ, IOBJ, MFIX, MBIAS) 

Called By : POTENT 

Purpose : Potentials for any conductors designated as fixed- 

potential or biased potential conductors are removed as un- 
knowns from the set of linear equations solved by the conju- 
gate gradient potential iteration system. Following converg- 
ence upon all other potentials of the system, proper biases 
are added to conductor I's new potential and the resultant 
biased potentials are placed in the appropriate places in 
arra^ PCOND (this is done in subroutine POTENT) . The only 
unknowns remaining in the system at this point are the new 
amour to of charge on the biased or fixed conductors. QCONCP 
finds the amount of charge on conductor i according to the 
original equation for this conductor: 



C. .4> + 

l] ] Fi F 


QCONCP calculates the charge on each conductor and 
places it in the appropriate location m array QCOND. 
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RUPDAT 

Calling Sequence : SUBROUTINE RUPDAT (ALPHA, IR, IAUN, ISPARE, 

BUF1, BUF2 , NG, NX, NY, NZ , RDOTR, RDOTRS , 
NC, RCOND, AUNCON, RDOTAU , MBIAS, MFIX) 

Called By ; POTENT 

Purpose : RUPDAT updates the' conjugate gradient vector r on 

file IR according to the equation 

r n+1 . r ” - a*(Au) n 

_ I 1 

The updated r data is written onto file ISPARE one grid at 
a time, and when all data has been updated, the unit numbers 
IR and ISPARE are switched. 

• x. 

RUPDAT also updates the conductor r values in array 
RCOND, and computes the two dot products r n+ ^*r n+1 and 
r n+1 -(Au) n . 
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SPACE 

Calling Sequence : SUBROUTINE SPACE (P, BUF3 , NX, NY, N2, 

W4 , UDOTAU) 

Called By ; ORCALL 

Local Variables : P: Storage space BUF2 

Purpose : Using the trilir.ear interpolation scheme, SPACE 

computes coefficient product values for grid points in the 
interior space (neither on inner edge n or outer edge) of 
the outer grid of the current pair. 

The interpolation requires both inner grid (PI) and 
outer grid (P2) data for u n (or p°) . Coefficient product 
values are placed in BUF3 . 
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UDOTUC 


Calling Sequence : SUBROUTINE UDOTUC (IU, UDOTU, BUF1, NG, 

NX, NY, NZ, NC, UCOND) 

Called By : POTENT 

Purpose: UDOTUC computes the dot product u n *u n for use in 

computing the IALG = 1 monotonic residual dot product 

,2 

r*r = u*u/c . 
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UDOTUP 

Calling Sequence ; SUBROUTINE UDOTU? (UDOTAU, TPROD ) 

Called By ; OPLANZ 

Local Variables ; TPROD: Coefficient product value for the 

grid point currently being treated 
by OPLANZ times u n for this point. 

Purpose : UDOTUP simply adds in the current term to u 11, (Au) n . 

In all other routines this simple operation is done internally, 
but in OPLANZ it caused compiler problems due to the over- 
flowing of an internal compiler table (UDOTAU is a double- 
precision variable) . 
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URSETO 

Calling Sequence : SUBROUTINE URSETO (ROUS, IU, IR, IAUN, 

BUF1 , BUF2 , NG, NX, NY, NZ , RDOTR, NC, 
AUNCON, UCOND, RCOND , ROUSCN, MFIX, MBIAS) 

Called By ; POTENT 

Purpose ; Given the set of linear equations for the potential 
of the system, Ap = ROUS, and given the coefficient product 
Ap° of A with the initial guess for the potentials, p°, URSETO 
defines the u° and r° arrays on files IU and IR according to 
the equation 

u° = r° = ROUS - Ap°. 

URSETO also calculates the initial dot product r°*r°. 
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UUP DAT 

Calling Sequence : SUBROUTINE UUPDAT (C, IU, IY, IR, IROUS, 

ISPARE, BUF1, BUF2 , BUF3, NG, NX, NY, NZ , 
NC, UCOND , YCOND , RCOND, ROUSCN, MBIAS, 
MFIX) 

Called By : POTENT 

Purpose : For algorithm 1, UUPDAT updates the conjugate 

gradient vector u on file IU according to the equation 

u n+1 = c*ROUS - y 

The updated u n+ ^ data is \*ritten onto file ISPARE one 
grid at a time and when all data has been updated, the unit 
numbers IU and ISPARE are switched. 

•UUPDAT also updates the conductor u values in array 

UCOND. 
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UUPDAU 

Calling Sequence ; SUBROUTINE UUPDAU (BETA, ISPARE, BUF1, 

BUF2 , NG, NX, NY, NZ, NC, UCOND, RCOND , 
MB IAS, MFIX) 

Called By : POTENT 

Purpose : For algorithm 0, UUPDAU updates the conjugate 

gradient vector u on file IU according to the equation 


u n+1 _ r n+ l + B , u r> 


The updated u n+ ^ data is written onto file ISPARE one 
grid at a time and when all data has been updated, the unit 
numbers IU and ISPARE are switched. 

UUPDAU also updates the conductor u values in array 

UCOND. 
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YINIT 

Calling Sequence : SUBROUTINE YINIT (IY, IAUN, BUF1, NG, NX, NY, 

NZ, NC, YCOND, AUNCON, MFIX, MBIAS) 

Called By : POTENT 

Purpose : Initialize Y-vector data on IY file and YCOND 

array for IALG = 1 potential algorithm. 
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YUPDAT 


Calling Sequence : SUBROUTINE YUPDAT (BETA, IY, IR, IROUS, 

ISPARE, BUF1 , BUF2, BUF3, NG, NX, NY, NZ , 
NC, YCOND , RCOND, ROUSCN, MB IAS, MFIX) 

Called By ; POTENT 

Purpose : For algorithm 1, YUPDAT updates the conjugate 

gradient vector y on file IY according to the equation 

y n+1 = g* y n + rouS - r n+1 

The updated y n+ ^ data is written onto file ISPARE 
one grid at a time, and when all data has been updated, the 
unit numbers IY and ISPARE are switched. 

YUPDAT also updates the conductor y values in array 

YCOND. 
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8. CODE DESCRIPTION - CHARGE SECTION 


The CHARGE section of NASCAP makes use of all the 
information — potentials, geometry, environment and 
material properties — developed elsewhere in the NASCAP 
code. Also, in this section is most of the physical 
content. 

CHARGE, the primary routine of this section, ob- 
tains the required information from external files, per- 
forms the LONGTIMESTEF option, and updates the charge 
density arrays. CALFLX calls routines to calculate the 
incident and emitted fluxes for each surface cell. Surface 
cell properties are decoded by GETCEL. For the test tank 
case, the incident flux is computed by SOURCE and ETNGUN; 
for the reverse trajectory case by PUSH. ELFLUX, PRFLUX, 
BKSCAT and SOLFLX compute emitted fluxes of electron- and 
proton-produced secondaries, backscattered electrons and 
photoelectrons, respectively. IOFLX performs printouts of 
these properties as requested. SHECAL (called from CHARGE) 
and PHOCUR perform the photosheath calculation. EFIELD 
and BFJELD provide the force fields required by the 
partic] e-tracking routines PUSH, ETNGUN and PHOCUR. 
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8.1 


PRIMARY ROUTINES - CHARGE, CALFLX 


CALF LX 

Purpose : Calls SOURCE, PUSH, ETNGUN , ELFLUX, BKSCAT , 

PRFLUX, GETCEL, IOFLX. (See flow chart. Figure 8.1.) 

Calling Sequence : SUBROUTINE CALFLX (ITYPE, PCOND, FLUXES, 

CNET, FLOW, IOBPLT, IPART) 


ITYPE - 
PCOND - 
FLUXES - 

CNET 

FLOW 

IOBPLT - 
IPART - 

Called By: CHARGE 


HASCAP operating mode 

Array of conductor potentials 

Array of net fluxes to surface 
cells 

Net charging current 

Low electron emission for each 
surface cell 

LUN of OBJDEF output file 

LUN of particle trajectory plot 
scratch file 


CHARGE 

Purpose : Primary routine for charge accumulation, transport 

and sharing for the NASCAP code. Performs LONGTIMESTEP 
option. (See flow chart. Figure 8.2.) 

Calling Sequence : SUBROUTINE CHARGE (NNX, NNY, NNZ, NIJG, 

ITYPE, QCOUD, PCOND) 

NNX, NNY, NNZ, NNG - Same as NX , NY , NZ , NG 

ITYPE - HASCAP operating mode 

QCOND - Array of charges on 

conductors 

PCOND - Array of potentials un 

conductors 


Called By : TRILIN 
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I'icjure 0.2. CHARGE flow chart 
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8.2 TEST TA11K ROUTINES 


ETNGUN 

Purpose : Calculate current density near object from tank 

beam. 

Calling Sequence : ETNGUN (IOBPLT) 

Input IOBPLT is the file which contains 
the object outline for plotting. 

Output COMMON/TANK/CURRENT (3,17,17), 
IBMCAL, COM, CURRENT (K, I, J) = J R (I,J) 
the Kth component of the current den- 
sity at grid location x = I, y = J. 

Called By : SOURCE 

Method : ETNGUN initially normalizes the net beam current to 

the input value. This is done by first integrating the cur- 
rent density and then scaling the input by the ratio of the 
desired current to integrated current. Next, the input 
spatial current density pattern is converted to current 
density as a function of angle at the beam source. For a 
given target size, only that part cf the beam necessary to 
completely cover the target is followed. The determination 
of the angular spread necessary for this is calculated by 
a Rutherford scattering type formula. Finally, test 
particle trajectories are followed forward in time until 
they either hit the target or exit the inner mesh. If they 
hit the target, their contribution is added to the current 
density array. 
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MAGCMP 

Purpose : Called by ETNGUN to correct velocity of electron 

for presence of constant magnetic field B during calibra- 
tion. 

Calling Sequence : SUBROUTINE MAGCMP (V, VO, RK, E, B) 

V{3) - uncorrected velocity, to be 
returned corrected 
VO - speed in code units 
RK - Z-distance (normal) from source 
to calibration plane 
E - electron energy (eV) 

B(3) - constant magnetic field 

Called By : ETNGUN 

SOURCE 

Purpose : Performs logical and I/O functions relevant to 

test tank particle tracking. 

C alling Sequence : SUBROUTINE SOURCE (IOBPLT) 

IOBPLT - LUtl of OBJDEF output file 

Called B y: CALFLX 
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8.3 


REVERSE TRAJECTORY ROUTINES 


FREAD 

Purpose : Reads DeForest's ATS-5 particle flux data. The 

unit number hour and day of the data are read off of IFLUX. 
JJJie hour and day must be precisely as appears on the graphics 
in the environmental specification portion of the final re- 
port. If this is not true, there will be an attempt to read 
past the end of file on IUNIT. Appropriate averages are 
taken of the data to determine mean energies and densities. 

Calling Sequence : SUBROUTINE FREAD (IFLUX) 

Called By : FLXDEF 

Local Variables : IUNIT: Unit where flux data is stored 

RHOUR: Time of requested data 

IDAY: Day of requested data 


FSPACE 

Purpose : FSPACE returns the phase space density from 

DeForest's data given the energy, velocity vector and 
species of a particle. 

Calling Sequence : FUNCTION FSPACE (ENG, V, I SPEC) 

Called By : PUSH 

Local Variables : C2: Absolute value of cosine of angle 

between velocity vector and standard 
vector 

S2: 1 - C2 
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PUSH 


Purpose ; For a given cell, calculate the incoming particle 
flux from a given ambient environment by the reverse tra- 
jectory technique. All test particles are emitted from the 
center of the cell. The selection of energies and angles 
are in common block FLUX. The actual flux densities are 
from FSPACE. The local flux is filled into F. EPOT(I) is 
the minimum final kinetic energy that a particle of type I 
can have upon impacting the object. That is, EPOT(I) = 0 
if the object repels type I particles and EPOT equals the 
absolute value of the local potential if the object attracts 
particles of type I. (See Figure 8.3.) 

Calling Sequence : SUBROUTINE PUSH (F, EPOT, ICELL, IPART) 

Called By : CALFLX 

Local Variables : X(3): Particle position 

X0{3): Center of the cell 

V^3)j Particle velocity 


SETAH 

Purpose ; Set up three standard 3x3 rotation matrices 
when the reference directions are the normals to <100>, 
<110> and <111> crystallographic planes. These are used 
in EFIELD and PUSH. 

Calling Sequence ; SUBROUTINE SETAN 
Called By ; PUSH 
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SETWE 


Purpose : For a given number of energy bins, NENG, calcu- 

lates the energies and weights such that each test particle 
would have equal weight if the incident flux were Maxwellian 
(E exp(-BE)). Also gives exact result for exponential flux 
(exp (-0E) ) . 

Calling Sequence : SUBROUTINE SETWE 

Called By: FLXDEF 


SMFDEF 

Purpose : Smooth raw ATS-6 data by simple linear blending 

to increase accuracy of average particle fluxes. The 
variable N sets how many data points on each side of a 
given data point will be included in the blended result. 

Calling Sequence : SUBROUTINE SMFDEF (FDEF) 

Called By : FLXDEF 

SETMAX 

Purpose : Constructs data for reverse trajectory sampling 

of a Maxwellian environment. 

Calling Sequence : SUBROUTINE SETMAX (A, FDEF) 

Called Bv: FLXDEF 
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8.4 


FLUX AND EMISSION ROUTINES 


BKSCAT 

Purpose ; Computes flux of backscattered electrons through 
Function BSCAT. 

Calling Sequence : SUBROUTINE BKSCAT (ICELL, ITYPE, FIN, FOUT) 

Called By; CALFLX 


BSCAT 

Purpose : Computes backscatter coefficient for electrons of 

energy (keV) incident normally (IOPT =1), at angle THETA 
(IOPT = 2) or isotropically (IOPT = 4). 

Calling Sequence : FUNCTION BSCAT (E, THETA, IOPT) 

E - Electron kinetic energy (keV) 

THETA - Incident angle (for IOPT = 2) 

Called By ? BKSCAT 
ELFLUX 

Purpose : Calculates incident electron flux and flux of 

electron-produced secondaries using subroutine ELSEC. 

Calling Sequence : SUBROUTINE ELFLUX (ICELL, ITYPE, FIN, FOUT) 

ICELL - Current surface cell 

FIN - Incident electron flux (output) 

FOUT - Flux of secondary electrons (output) 
ITYPE - NASCAP operating mode 


Called By : 


CALFLX 
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ELSEC 

Purpose : Calculates secondary yield for electrons of 

energy E (keV) which are normally incident (I0PT=1) 
incident at angle THETA (I0PT=2) , or isotropically inci- 
dent (I0PT=4). Subroutine ELSEC1 is identical to ELSEC. 

Calling Sequenc e: SUBROUTINE ELSEC (YIELD, E, THETA, IOPT) 

Called By : ELFLUX 

IOFLX 

Purpose : Output of fluxes, potentials and fields across 

dielectrics for selected surface cells. 

Calling Sequence : SUBROUTINE IOFLX (EIN, EOUT, PIN, POUT, 

EBACK, PCUR, PCND, ICELL) 

Called By: CALFLX 


PRFLUX 

Purpose : Calculates flux of incident protons and resulting 

secondaries using subroutine PROSEC. 

Calling Sequence : SUBROUTINE PRFLUX (ICELL, ITYPE, FIN, FCUT) 

Called By : CALFLX 

PRO SEC 

Purpose : Calculates secondary electron yield for incident 

protons of energy E (keV) , incidence nay be normal (I0PT=1) , 
isotropic (IGPT=4) or at angle THETA (I0PT=2) . 

Calling Seque nce: SUBROUTINE PROSEC (YIELD, E, THETA, IOPT) 

Called By : PPELUX 
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SOLFLX 

Purpose : Calculates flux .of photoelectrons from surface cell 

using material properties and shadowing calculation output. 

Calling Sequence : SUBROUTINE SOLFLX (ICELL, PCUR) 

Called By: CALFLX 


8.5 FORCE CALCULATION, VECTOR MANIPULATION AND INFORMATION 
DECODING ROUTINES 


ADDVXB 

Purpose : Computes vector product 

VOUT = VIN + CMAG*VIN * B 

Calling Sequence : SUBROUTINE ADDVXB (VOUT, VIN, B, CI1AG) 

Called By : . PUSH, PHOCUR, ETNGUN 

BFIELD 

Purpose : Computes magnetic field (in W/m^) B at point X. 

Calling Sequence ; SUBROUTINE BFIELD (3, X) 

Called By : PUSH, PHOCUR, ETNGUN, MAGCMP 
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EFIELD 


Purpose : Compute the local electric field components for a 

plasma particle and determine whether the particle has 
collided with the object or exited the outermost grid. 

Calling Sequence : EFIELD (E, X, ITEST) 

Input X(3) is the particle position. 

Output E(3) is the vector electric field 

1 0 particle in free space 

1 particle is outside the second mesh 
-1 particle inside the object 

Called By ; PUSH, PHOCUR and ETNGUN. (A flowchart of EFIELD 

is included in the description of subroutine PUSH, 
Figure 8.3.) 

Method : EFIELD calls the following subroutines: 

LCODE - decodes the element table for cells that 
are partially full. 

TRAND - transforms the fractional particle position 
inside a special cell to the- standard ori- 
entation of that cell (e.g., from 101 to 
110 for weages) . 

INSIDl - determines (using the transformed particle 
position) whether the particle is inside 
or outside the object. 

TRANA - transforms the matrix of eight potentials 
of a special cell to the standard orienta- 
tion of that type of cell. 

ESPEC - actually calculates all the electric fields. 
This is dene in the appropriate standard 
orientation for special cells. 

TRANE - transforms the special cell electric field 
vector from the standard orientation back 
to the orientation of the particular element 
that the particle is in. 


234 



I . ■ — l 


I 


ESPEC 

Purpose : Calculates local electric field in a zone in the 

standard orientation. (See Figures 6.19 to 6.24.) 


Calling Sequence : ESPEC (E, X, P, ITYP) 


Input X(3) position in cell. 

P(8) potential at the eight vertices of the cell 
(for some cells not all the potentials are 
used) . 

= 0 cube (empty cell) 

= 1 wedge (110) 


ITYP 

(cell type) 


= 2 cube with line on z face 
= 3 tetrahedron (111) 

= 4 truncated cube (111) 


Output E (3 ) local electric field 


Called By ; 

Method : Uses analytical derivative of the trilinear (cube) , 

bilinear (wedge) or linear (tetrahedron) basis sets. The 
formulation for the truncated cube is only approximate and 
uses trilinear weights with an extrapolated potential at the 
missing vertex. 


FILINP 

Purpose : Fills in the boundary region of the outer grid with 

values from the inner grid. 

Calling Sequence : SUBROUTINE FILINP (PI, P2, NX, NY, NZ, 1JG) 

Called By : TRILIN, CHARGE 
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GETCEL 

Purpose : Places geometrical and electrostatic properties of 

surface cell in common block SRFPRM. 

Calling Sequence : SUBROUTINE GETCEL (ICELL, PCOND) 

Called By : CALFLX, SHECAL 


INSID1 

Purpose : Determines whether transformed point DX is interior 

given cell of type ITYP in standard orientation. 

Calling Sequence : FUNCTION INSIDi (DX, ITYP) 

Called By : EFIELD 


LCODE 


Purpose : Decodes entry in LTBL (volume cell element table) . 

(See EFIELD.) 

Calling Sequence : SUBROUTINE LCODE 

L - Element table entry 

ITYP - Cell type 

I0R(3) - Orientation code 

JCURCL - Index entered by NUMLTB 


Called By : EFIELD 


TRANA 

Purpose : Rearranges the eight vertices of a cube AI (2,2,2) 

into AJ(2,2,2) by the same transformation of axes that trans- 
form vector 1(3) into J(3). 

Calling Sequence : SUBROUTINE TRANA (AI, AJ, I, J) 

Called By: EFIELD 
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TRAND 


Purpose ; Transforms real positive differential displacement 
vector DI (3) into DJ(3) as vector 1(3) is transformed into 
JO). All components of DJ are positive and less than unity. 

Calling Sequence : SUBROUTINE TRAND (DI, DJ, I, J) 

Called By: EFIELD 


TRANE 

Purpose : Transforms electric field vector El (3) into EJ(3) 

as vector 1(3) goes to J (3) . 

Calling Sequence : SUBROUTINE TRANE (El, EJ, I, J) 

Called By : EFIELD 


TRAAN 

Purpose : Takes a rotation matrix AI and transforms it into 

AJ by the same axis interchanges that converted vector 1(3) 
to J(3). It determines the interchanges as well as transforms 
the matrix. 

fa lling Sequence : SUBROUTINE TRAAN (AI, AJ, I, J) 

Called By : PUSH 


VADD 

■*> 

Purpose : Add two 3 vectors: A = B *■ C. 

Calling Sequence : SUBROUTINE VADD (A, B, C) 

Called By : Various routines 
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VADDS 

Purpose : Add two 3 vectors after multiplying one by a scalar: 

A = 3 + SCALAR x C. 

Calling Sequence : SUBROUTINE VMULT (Y, A, X) 

Called By: Various routines 


VMULT 

Purpose : Multiply a three - vector , X, by a 3 x 3 matrix A, 

and yield result Y: Y = AX. 

s* 

Calling Sequence : SUBROUTINE VMULT (Y, A, X) 

Called By: Various routines 


VSET 

Purpose : Set 3 vector A equal to 3 vector B: A = B. 

Calling Sequence : SUBROUTINE VSET (A, B) 

Called By : Various routines 

8.6 PHOTOSHEATH ROUTINES 


Ph OCUR 

Purpose : Photoelectron and secondary electron contributions 

to space charge densities and to surface currents for a given 
zone are calculated in this routine by following test particles 
forward in time. The particles are pushed in a fashion entirely 
analogous to PUSH. If particles are within one zone of a sur- 
face cell, all three components of the current are stored. 

This is used to define local particle fluxes. 
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Calling Sequence : SUBROUTINE PHOCUR (PHOJ, QEMIT) 

Called By : SHECAL 

Local Variables : XO(3): Center of cell 

XI (3,4): Places where test particles are 

emitted 

SHECAL 

Purpose : Performs logical functions and charge sharing for 

photosheath calculations. (See flow chart. Figure 8.4.) 

Calling Sequence : SUBROUTINE SHECAL (PHOJ, FLOW, FLUXES, 

PCOND) 

PHOJ (3, 1024) - Array for fluxes in 

selected volume cells 

FLOW - Low energy electron flux 

emitted from each surface 
cell 

FLUXES - Net flux to each surface 

cell 

PCOND - Conductor potential array 

Called By: CHARGE 
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Figure 8.4 


Photosheath calculation (SHECAL) flow chart 
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9. CODE DESCRIPTION - SATPLT, HIDCEL AND 
OTHER GRAPHICS SECTIONS 

SATPLT calls all routines responsible for satellite 
illustration plots. HIDCEL creates three-dimensional pers- 
pective plots, while MATPLT produces the chaded area sil- 
houettes. CONPLT plots potential contours, ROUSP plots ROUS 
charge contours, and PARPLT plots particle trajectories. 

Only the three-dimensional hidden-line plots are com- 
plex enough to warrant an algorithm description. Following 
is a description of the HIDCEL algorithm and documentation of 
all plot-related subroutines. 
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9.1 


SHADOWING ALGORITHM DESCRIPTION 


For purposes such as calculating photoemission effects 
it is desirable to be able to calculate the shadowing effects 
of the satellite upon itself for an arbitrary angle of inci- 
dent sunlight. A subroutine called IIIDCEL computes for each 
surface grid cell the area of that part of the cell which re- 
mains exposed when partially or completely shadowed by some 
other portion of the satellite (or, of course, when not 
shadowed at all). The program is completely general: For 

any satellite which NASCAP can treat, any incident angle may 
be specified as a vector from the center of the satellite 
out towards the viewer (sun) . In addition to calculating the 
cell-by-cell areas needed for code computations, HIDCEL plots 
a diagnostic (and informative) hidden line .3-D plot of the 
satellite as viewed from the incident angle requested. 

The code accomplishes both the area calculations and 
3-D plots by a long series of polygon maskings. The set of 
cell faces whose normals produce positive dot products with 
the vector pointing towards the viewer are selected and pro- 
jected into 2-D polygons. We call this set of individual cell 
polygons the Al polygons. A second but much shorter set of 
polygons is selected and projected ;nto 2-D by finding those 
faces of each parallelepiped, wedge, tetrahedron, octagon, 
etc. forming the satellite (as determined by user input), 
which meet the dot product criterion. We call this set the 
A2 polygons. 

The cell-by-cell (Al) set of 2-D polygons is then 
masked one at a tine against each of the larger "building- 
block" (A2) polygons, and at each step we discard any poly- 
gon which is partially or completely shielded and replace it, 
if necessary, with the polygon (or polygons) representing 
that part of the original cell which remains visible. 
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An A1 polygon which is completely shadowed, but is 
masked by several different surfaces, will then be "whittled 
down" a bit at a tine, each masking leaving an appropriately 
smaller polygon until the last masking for that cell finally 
removes it from the system. A polygon whose niddle is 
masked but whose ends remain visible will be split into two 
or more new polygons; the code maintains records, however, 
indicating which original cell has generated each new polygon. 
Consequently, vrhen all naskings have been completed and a 
finalized set of nasked A1 polygons exists, we may calculate 
areas of each polygon and add into the appropriate cell area 
the area of each polygon contributing to that cell. 

The areas stored for further NASCAP calculations are 
fractional areas with the dot product (2-D projection) effect 
removed. We do this by calculating the areas of each cell's 
2-D projection before any masking takes place and then comput- 
ing the fraction of this area which remains visible following 
all naskings. This procedure also precludes any area 
normalization problems. 

Figures 9.1 through 9.3 illustrate plot output of 
the shadowing program for several different angles of viewing. 
Figures 9.4 and 9.5 demonstrate printed output of the areas 
and fractional areas for one case. 

As is evident, for a satellite of this complexity, we 
may be masking on the order of 500 cell faces against on the 
order of 150 building-block faces for a total of about 75,000 
individual markings. Computing time (including plot genera- 

3 

tion) on S 's UNIVAC 1108 varies from 5 to 60 seconds per 
view, depending on satellite complexity. 
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Plot output from shadowing routines. 
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Figure 9.3. Plot output from shadowing routines. 
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Figure 9.4. Printed output from shadowing routines. 
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9.5. Printed output from shadowing routines. 
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F Kjuie 9.3. Continued 



9.2 


SUBROUTINE DESCRIPTIONS 


ADDAI 

Calling Sequence - SUBROUTINE ADDAI (A3, N3, Al, NAl, NV1, 

IFFS, IFF1 , AIRS, AIR, A1BS, A1FS, A1B, 
A1F, XMIN1, XMAXl, YMIN1, YMAX1 , IA1) 


Called By : HIDCEL 

Purpose ; ADDAI adds a new polygon and its associated in- 
formation to the end of the Al polygon list. 
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AREA 

Calling Sequence : FUNCTION AREA (A, N) 

Called By : HIDCEL, F INVER 

Local Variables : A(2,N): Vertices of a 2-D polygon 

N: Number of vertices 

Purpose ? AREA'S value is the area of polygon A. 


251 



A1C0MP 


Calling Sequence : SUBROUTINE A1COMP (Al, NA1 , HV1, IFF1, AIR, - 

A1B, A1F, XMIN1, XMAXl, YMIN1, YMAXl, IREMOV) 

I 

Called By : KIDCEL 

Purpose : A1COMP is called following masking of all Al poly- 

gons against one A2 polygon. A1C0MP removes all flagged Al 
polygons (any which have been partially or wholly shielded) , 
and compresses the list by inserting polygons from the end 
of the list into those spaces occupied by the removed ones. 


252 



I 


t 


t 


I 


A1GEN 

Calling Sequence ; SUBROUTINE A1GEN (DIR, Al, A1B, A1F, 

NA1, NV1, DV1 , AIR, IFF1, XM1N1 , XMAXl, 
YMIN1, YMAXl, NSURF, JSURF, IP) 

Called By ; HIDCEL 

Purpose : A1GEN generates the vertices in 3-D of all Al 

individual cell faces whose normals have positive dot 
product with the direction of viewing. Each of these 3-D 
polygons is then projected into 2-D by calling PRSPPJ. The 
final set of 2-D vertices is stored in the Al array. In 
addition, 2 -d XllIN's, XWAX's, YMIN's, and YMAX's are cal- 
culated and stored for use in rough masking and behindness 
tests. 
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A2GEN 

Calling Sequence : SUBROUTINE A2GEN (DIR, A2, A2B, A2F, 

IFAC , NA2 , IFF2 , DV2, A2R, NV2, NBLK, 
NPBLK, IPER, XMIN2, XMAX2, YMIN2 , 

YMAX2, IP) 

Called By : HIDCEL 

Purpose : A2GEN generates che 3-D and 2-D projected co- 

ordinates for the "building-block" face A2 surfaces whose 
normals have a positive dot product with the direction of 
viewing. A2GEN also computes and saves minima and maxima 
for use in behindness tests. 

Internal Organization : We start with a loop over all 

building-blocks making up the satellite as defined by user 
input to OBJDEF . The number of points in the block 
(NPBLK (IBLK) ) tells us which type block it is of the five 
possible types: rectangular, parallelepiped, wedge, 

octagon, or "FILlll" type. From the type, we know how 
many faces this block has, and we loop over each face- 
the current block. 

Using face normal data stored in the IFAC array 
(see SETFAC) , we find the normal for the current face and 
place it in IFF (3). We compute the dot product of this 
face normal with the direction of viewing, and skip this 
face if its dot product is less than or equal to l.E-10. 

If we keep this face we must f ine*-Mk£»- -3-D vertices. 
We pick tne vertices for this face in order out of array 
IPER according to the face index, the type of block, and 
the vertex index data stored in IV (see SETFAC) . 

Next, we project the 3-D vertex set into 2-D by 
calling PRSPPJ for each vertex. We also store the distance 
to the viewer for each vertex in DV2. The closest and 
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furthest points to the viewer are found and stored in A2F 
and A2B. The 2-D minima and maxima are calculated and 
stored, and the first vertex is appended to the vertex 
list. 

All faces of all blocks are looped through in this 

manner . 


\ 
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A2PL0T 


Calling Sequence : SUBROUTINE A2PLOT (N3LK, NPBLK, IPER, A2P) 

Called By : HIDCEL 

Purpose ; A2PL0T plots the non-hidden line plots of the 
building-block faces when HIDCEL is called with IAREA = 1. 

If drastic errors have been made in defining the building 
blocks which make up the satellite, A2PL0T plots will re- 
veal then in an obvious manner. 
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CONPLT 

Calling Sequence : SUBROUTINE CONPLT (IXYZ, IUD, NX, NY, NZ, 

PI, P2, IPER, PER, IOBPLT ) 

Called By : TRILIN 

Local Variables ;. IXYZ: IXYZ =1, 2 or 3 specifies contours 

through a fixed x, y or z plane 

IND: Fixed value of x, y or z; index in 

next-to-inner grid at which cut is 
to be taken 

PI: Inner grid potentials 

P2: Next-to-inner grid potentials 

Purpose : CONPLT plots potential contours for the inner two 

grids. The 3-D plane at which the potentials are plotted is 
specified by a fixed x, y or z. IXYZ determines which dimen- 
sion is fixed, and IND is the fixed value. 

IND is an index in the outer grid of the pair. 

Internal Organization : CONPLT fills in a two-dimensional array 

P with the potential values at the specified cut. P has the 
fineness of the inner grid, but encompasses both grids. Poten- 
tial values at the new finer points in the outer grid area of 
P are found by linear interpolation, while the values for 
the inner grid are simply those already present in PI. 

CONTOR is called to plot the potential contours and 

titles. 

The satellite silhouette is plotted by reading in 
the building-block data from file IPLT = IAES (IOBPLT) and 
plotting the block perimeters in IPER = PER projected into 
the appropriate plane, and translated according to the new co- 
cordinate system. 
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Calling Sequence : 


SUBROUTINE CONPOT (BUF1, NX, NY, NZ, PCOND) 


Called By ; TRILIN 

Purpose ; CONPOT sets the potentials at points in the interior 
of the satellite (otherwise untreated) equal to the potential 
on conductor one. 
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Calling Sequence : SUBROUTINE F INVER (IHEW, NEW, A3, A2, N2, 

N3, 113INT, NA2SG3, IA2SG3, INS2, IN2, NJOIN, 
I JOIN, IP) 

Called By : SHIELD 

Local Variables : AeWs Number of points currently 

in the new polygon 

NEW: Index of new polygon being 

considered 


N2: Humber of vertices in the 

current A2 polygon 


H3(12): Number of vertices (final) 

in each of the new polygons 

N3INT(12): Number of intersection points 

for each of the new polygons 

NA2SG3 (2,12) : Number of A2 segments inter- 

sected by up to two inter- 
section points per rew poly- 
gon, for up to twelve new 
polygons 


IA2SG3 (2, 2, 12 ) : Indices of A2 segments inter- 

sected for the intersection 
points indexed by NA2SG3 


INS 2 (9) : INS2 (I) is 0 or 1 if the Ith 

vertex of A2 is outside or 
inside A1 

IN2: Number of A2 vertices inside 

A1 


NJOIN: Not used 


I JOIN: Not used 

IP: Print flag 


Purpose : FINVER completes the vertex set of a new polygon 

consisting of 1 or more A1 vertices and two points where A1 
intersects the A2 perimeter, by adding those A2 vertices which 
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fall between the two intersection points and which are inside 
the original A1 polygon. 

Interna? Orga ni zation ; First we determine if the first and 
last points of the new polygon are the same. If they are 
not, we proceed to the normal loops of the routine. If 
they are, we check if INEW =3. If so, we have tv/o points 
different to the finer accuracy of the intersection routine, 
hut the same to the coarser accuracy of PTCOMP. We keep 
the finer accuracy polygon and return. If INEW / 3, we have 
a redundant last vertex and remove it to yield an already 
complete polygon. We return with this polygon. 

If the two intersection points are on the sane A2 seg- 
ment, we are done already. 

The normal section of the routine marches around the 
vertices of A2 between the two intersection points, first in 
one direction, and then, if necessary, in the other. It 
searches until it finds a complete set of A2 vertices between 
the two points which are all within Al. This is the set of 
points we need. If any point in the first direction is out- 
side Al, we skip out and search in the other direction. If 
no such set is found in the second direction, we error return. 

For some very complex cases, we may need to calculate 
the area of the new polygon taken with each of the two pos- 
sible vertex sets. The set providing the lesser area is 
the correct set. 
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HIDCEL 


Calling Sequence : SUBROUTINE HIDCEL (PDIR, IPLOT, NX, NY, 

NZ, I AREA) 

Called By : NASCAP , SATPLT 

Local Variables : PDIR(3) : Vector (not necessarily normalized) 

indicating direction from center of 
satellite to viewer {or sun, for 
shadowing) . 

Purpose : HIDCEL calls all routines associated with the 3-D 

plots. There are two modes of calling HIDCEL: (1) IAREA = 0 

implies that this call is merely for satellite illustration, 
and both a "building-block" non-hidden line plot and a cell- 
by-cell hidden line plot will be plotted. In this case, no 
areas will be calculated; (2) IAREA jt 0 implies that this call 
is to generate hidden-cell fractional areas for each of the 
(up to 1024) surface cells. These areas will be written onto 
unit IAREA at the end of HIDCEL. A cell-by-cell hidden line 
plot is also plotted in this case, but no building-block plot 
is generated. 

Internal Organization : HIDCEL first normalizes the PDIR (3) 

vector pointing from the center of the satellite towards 
the viewer. The normalized data is placed in DIR (3). PPJSET 
is then called to set up the 3-D to 2-D (PRSFPJ calls) and 
2-D to 3-D (PRSPPK calls) conversion matrices. Tata is read 
in from file IPLOT = IABS(IOBPLT) holding the building-block 
information which defines the satellite. If we are in node 
IAREA = 0, subroutine A2PL0T is called to generate the 3-D 
non-hidden-line plot of the building blocks. 

Next, A2GEN is called to generate all A2 polygons and 
all data associated with them. The set of A2 polygons is the 
set of building-block faces whose normals have a positive dot 
product with DIR, projected into the 2-D focal plane perpendic- 
ular to the direction of viewing as set up by PPJSET. 
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A1GEN is called to generate the set of Al polygons 
and their associated data. The set of A1 polygons is the 
set of surface cell faces whose normals have a positive dot 
product with DIR, projected similarly into 2-D. The 2-D 
area of each A1 polygon is calculated and the information 
written onto file IAREA. These areas reflect the angle of 
a cell relative to the direction of view, and represent the 
entire 2-D area of each cell before shielding. 

Two nested loops march through every combination of 
each A1 polygon with each A2 polygon. The outer loop loops 
through the A2 polygons. Thus, for one building-block face 
polygon, each cell face polygon is checked against it one 
at a time. 

If the A2 vertex closest to the viewer is at least as 
far as the farthest A1 vertex, A2 is behind A1 and this com- 
bination is skipped. 

2-D XMIN's, XMAX's, YMIN’s, and YMAX*s are generated 
in AlGEM and A2GEN for each A1 and A2 polygon. If A2 ' s 
XMAX is greacer than or equal to Al's XMIN, or if A2's XI1IN 
is less than or equal to Al's XMAX, or similarly for the 
YIlIN's and YIlAX’s, we skip this combination. 

A variable called I21ASK, initially set to 0, indicates 
whether we have yet determined if Al is behind A2 and must 
be masked. If the A2 vertex furthest from the viewer is at 
least as close as the Al vertex closest to the viewer, Al 
is behind A2 and IMASK is set to 1. If IMASK is 0 at this 
point, however, several possibilities still exist as to 
whether the complex shielding algoritluu nusc be executed for 
this combination (e.g., Al may be on A2 in which case it 
must not be masked against A2). 'Subroutine SHIELD is called 
to make the final determination of whether masking is to be 
done, and to do it if so, or, if IMASK is already 1, to 
simply execute the shielding algorithm, skipping any further 
tests on behindness. 
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ISHLD is the main output indicator from SHIELD. 

ISHLD = -1 implies that no changes need be made in any of 
the polygon bookkeeping arrays. This occurs when it is 
determined that the polygons need not be masked against one 
another due to Al's being in front of or on A2, or when the 
polygons are masked against one another but there is no 
overlap. ISHLD >_ 0 implies that the polygons were masked, 
that there is overlap, and that ISHLD new A1 polygons have 
been created and are currently in array A3 waiting to be 
transferred to array Al. Thus, if ISHLD = 0, the current 
A1 polygon was completely shielded by A2 and zero new poly- 
gons have been created. Regardless of whether any new poly- 
gons have been created or not, the current Al polygon is 
flagged for removal by subroutine A1COMP by setting 
IREMOV (IAl) equal to 1. 

If new polygons have been created (ISHLD >_ 1), sub- 
routine ADDA1 is called for each of them to add their 
respective information to the tail end of the Al polygon 
arrays. 

Next we loop back around to the beginning of the Al 
polygon loop to begin checking of the next Al polygon against 
the current A2 polygon. 

When all Al polygons (except those new ones created 
by comparison to the current A2 polygon) have been checked 
against the current A2 polygon, subroutine AlCOriP compresses 
, . -w-sgjie arrays containing Al polygon information by removing each 
polygon flagged in IREIOV, and replacing its data with data 
from the end of the list (starting with the last new Al 
polygon created) . Thus, when checking against the next A2 
polygon begins, the list of Al polygons which we check it 
against is exactly what we want it to be; the most current 
set of Al polygons containing some already partially masked 
by previous A2's and containing no superfluous polygons. 
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When all A2 polygons have been checked against all A1 
polygons, we are left with a final set of non-overlapping A1 
polygons representing those parts of all surface cells re- 
maining visible when viewed from direction DIR. We loop 
over this finalized set of A1 polygons calculating an overall 
2-D XMIN, XliAX, YMIN and Y11AX. We also calculate the area 
of each polygon and add it into the appropriate location in 
array AREAl. The index in AREA1 is the index of the original 
surface cell whose 2-D projection^was the parent polygon of 
the current A1 polygon. This index is stored in bits 6-15 
of IFFl(IAl). A given surface cell can generate two or more 
partially shielded new polygons if some thin part of the 
satellite crosses over the middle of the surface cell. For 
diagnostic purposes, a message to this effect is printed out 
along with the Al polygon area whenever such a situation is 
encountered. 

Next, the original unshielded areas of all surface 
cells are read back in from unit IARCA. (Cells which did 
not face the viewer have this area already set to zero and 
hence, in this sense, have already been shielded.) A second 
loop over the surface cells computes the fractional areas 
for each cell: final shielded 2-D area/initial unshielded 

2-D area. In this way, the effect of the dot product on the 
area is removed. 

Next, a user area for plotting is set up using the 
2-D overall minima and maxima just calculated. A last loop 
over all final Al polygons plots the perimeter of each poly- 
gon by simply connecting its vertices. 

Finally, the IAREA file is rewound, the set of 1024 
fractional areas is written onto it, and IAREA is rewound 
again for use by the rest of the code. 
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Calling Sequence: FUNCTION 

Called By: SHIELD 

INSIDD (A, AD, P 

, PD, N, IP) 

Local Variables: A(2,10): 

Single-precision 

polygon vertices 

AD (2,10) : 

Double-precision 

polygon vertices 

P: 

Single-precision 

point 

PD: 

Double-precision 

point 

N: 

Number of vertices in polygon 

IP: 

Print flag 



Purpose : INSIDD determines if point P is inside polygon A. 

If P is on the perimeter of A, it is considered inside. 
INSIDD « 0 implies P is not inside A. INSIDD = 1 implies P 
is inside A. 


INSIDD is used instead of INSIDE for determining if 
vertices of a polygon are inside the other and for resetting 
them to be exactly on the polygon perimeter if they are off 
by only some epsilon. This resetting is accomplished by 
calling PROLIN to do a double-precision calculation finding 
the projection of point P on the segment of A which it is 
almost on. 
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INSIDE 

Calling Sequence ; FUNCTION INSIDE (A, P, N, IP) 

Called By : SHIELD 

Local Variables : A(2,10): Vertices of a polygon 

P (2) : Point 

N: Number of vertices of polygon 

IP: Print flag 

Purpose : INSIDE determines if point P is inside polygon A. 

If P is on the perimeter of A, it is considered inside. 
INSIDE = C implies P is not inside A. INSIDE = 1 implies 
P is inside A. 

INSIDE is a single-precision function. 
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INTCHK 

Calling Sequence : SUBROUTINE I1ITCHK (II, Al, AID, A2, A2D, 

INT, N2 , IHD1NT, X1NT, NA2SEG, IA2SEG, 

INS2, IP) 

Called By ; SHIELD 

Local Variables : II: Index of current Al polygon 

segment 

AID: Double-precision copy of Al 

A2D: Double-precision copy of A2 

INT: (output) number of inter- 

sections of segment II of Al 
with A2 perimeter 

N2: Number of A2 segments 

INDINT: Not used 

XINT (2,10) : (output) intersection points 

found 

NA2SEG (10 ) : (output) number of A2 segments 

intersected by up to 10 inter- 
section points 

IA2SEG (2,10) : (output) indices of A2 segments 

intersected by up to 10 inter- 
section points 

INS2(9): INS2(I) = 0 or 1 if vertex I 

of A2 is outside or inside Al 

IP: Print flag 

Purpose : Given the Ilth segment of Al, 7NTCHK checks for any 

and all intersections of this segment with all A2 segments. 
INTCHK also keeps track of which A2 segments are intersected, 
eliminates duplicate intersection points, and orders those 
kept according to how close to the start of the Ilth Al 
segment they are. 
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Internal Organization : I1JTCHK starts with a loop over the A2 

segments. INTSEC is called to find any intersections of the 
A1 segment with the current A2 segment. If there were none, 
we go on to the next A2 segment. If there was one, we incre- 
ment IUT and record the A2 segment information in IJA2SEG 
and IA2SEG (INTSEC has already placed the intersection point 
in XINT) . The distance of this point to the start of the A1 
segment is computed and stored in D for sorting purposes. If 
there has been at least one previous intersection, we check 
to see if we have a redundant intersection and eliminate it 
if so. 

If the A1 segment does not "intersect" the A2 segment, 
but partially or wholly overlaps it, we do not need to count 
the intersection but we must record the A2 segment information 
in NA2SEG and IA2SEG for use by F INVER. 

When all A2 segments have been checked, we sort all 
intersection points (if more than one) according to their 
distance from the start of Al, so that the closest ones will 
be treated first. 
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INTSEC 


Calling Sequence ; SUBROUTINE INTSEC (XI, X2, PXI, INTFLG, IP) 
Called By ; INTCHK 

Local Variables : XI (2,2): First double-precision line segment 

X2(2,2): Second double-precision line seg- 

ment 


PXI (2): (output) single-precision inter- 

section of XI and X2, if any 

INTFLG: (output) INTFLG = 0 implies no 

intersection; INTFLG = 1 implies 
an intersection; INTFLG = -1 im- 
plies the segments overlap. 

IP: Print flag 


Purpose : INTSEC determines if segment XI intersects segment 

X2, and if so, calculates the intersection in double-precision. 

INTSEC finds the equations of the two lines, finds 
their intersection, and determines if that intersection is on 
both of the line segments by calling PTLINE. 
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Calling Sequence ; SUBROUTINE LINPLN (XLl, DIR, XP, INB, XINT) 
Called By : SHIELD 

Local Variables ; XLl (3): Start point of line 

DIR(3): XLl + DIR forms finish point of 

line 


XP(3): Point on the plane 

'INB: Word containing bit-encoded normal 

to plane 

Bits 0,1: x-component of normal 

Bits 2,3: y-component of normal 

Bits 4,5: 2-component of normal 

XINT(3): (output) intersection of the line 

with the plane 

Purpose : LINPLN finds the intersection XINT of a 3-D line 

with a 3-D plane. 



I1ATCAL 


Calling Sequence : SUBROUTINE MATCAL (NX, NY, N3, IOBPLT) 

Called By ; SATPLT 

Purpose : MATCAL reads the satellite surface cell and material 

information from unit IOEPLT and calli. MATPLT to plot the 
satellite material composition silhouettes. 
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MATPLT 


Calling Sequence : SUBROUTINE MATPLT (NX, NY, NZ) 

Called By : MATCAL 

Purpose : MATPLT plots satellite silhouettes illustrating sur- 

face material composition via fifteen different shadings. The 
default node is to plot six views fron the plus or minus x, y 
and z directions, plotting all cells whose normal has a com- 
ponent in the direction of viewing. In cases where a cell 
face is partially or completely shielded by other cells in 
front of it, only the visible part of that cell (if any) is 
plotted and shaded. The sacrifice here is that inner sur- 
faces hidden by outer surfaces may never be seen via the de- 
fault vievrs. 

Hence, additional input allows the user to specify for 
any of the six basic viewing directions, additional plots show- 
ing only those cells within a specific range of y or z 
values according to the view involved. 
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M0VEA1 


Calling Sequence ; SUBROUTINE ilOVEAl (II, 12, Al, NV1, IFF1, 

AIR, A1B, A1F, XfllNl, XMAX1, YMIN1, YMAXl) 

Called By : A1COMP 

Purpose ; MOVEAl noves the information for the Ilth Al polygon 
into the I2th position in all the Al information arrays. 
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Calling Sequence ; SUBROUTINE PARPLT (IXYZ, I PART, IOBPLT, NX, 

NY, HZ) 

Called By ; TRILIN 

Local Variables : IXYZ: IXYZ = 1, 2 or 3 specifies particle 

trajectories projected onto the y-z, 
x-z or x-y plane (fixed x, fixed y 
- or fixed z) 

Purpose : PARPLT plots particle trajectories generated by sub- 

routine PUSH projected onto either of the y-z, x-z or x-y 
planes. Separate plots are made for electrons and protons 
(species one and two) . 

The particle data is read from file IPART and plotted 
by connecting each successive location for one particle to 
the next. 

The satellite silhouette is plotted by reading in the 
building-block data from file IPLT = IABS (IOBPLT) , and plot- 
ting the block perimeters in IPER = PER projected onto the 
appropriate plane. 
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PROLIN 

Calling Sequence : SUBROUTINE PROLIN (A, X, XP, IP) 

Called By ; INSIDD 

Local Variables : A (2, 2): Double-precision line segment 

X(2): Double-precision point 

XP(2): (output) projection of X on A 

IP: Print flag 

Purpose : PROLIN computes the projection XP of point X on 

line segment A. 
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PRSPPJ 


Calling Sequence : Entry Point PPJSET: SUBROUTINE PPJSET 

(X1Q, X2Q, X3Q, XIV, X2V, X3V, D) 

Entry Point PRSPPJ: - SUBROUTINE PRSPPJ 

(x, y, z , ::p, yp, zp) 

Entry Point PRSPPK: SUBROUTINE PRSPPK 

(XJ, YJ, ZJ, XI, YI, ZI) 

Called By : HIDCEL (PPJSET) 

A1GEN, A2GEN (PRSPPJ) 

SHIELD (PRSPPJ, PRSPPK) 


Local Variables : X1Q: 

X2Q: 


X3Q: 


XIV: 
X2V : 
X3V : 
D: 


x coordinate of satellite center 
(inner grid center) 

y coordinate of satellite center 
(inner grid center) 

2 coordinate of satellite center 
(inner grid center) 

x coordinate of canera lens 

y coordinate of canera lens 

2 coordinate of camera lens 

camera focal length 


Y . i coordinates of 3-D point to be 
1 projected into 2-D 
Z: » 

XP • \ 

I (output) coordinates of 2-D pro^ec- 
YP: , tion and (ZP) distance to viewer 

ZP: ' 


yj. ( (output) coordinates of 3-D reverse- 
’ I projection Doint 
ZJ: ) 

XX * 

J coordinates of 2-D point to be reverse- 
YI: . projected and (ZP) distance to viewer 
J of 3-D point to be found 
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Purpose : PPJSET sets up the 3-D to 2-D and 2-D to 3-D pro- 

jection matrices. PP.SPPJ projects from 3-D to 2-D; PRSPPK 
projects from 2-D to 3-D. The method is as follows: 

Imagine a camera properly oriented to take a picture of 

a collection of points and line segments. Let P v be the point 

at the center of the lens of the camera and let be a point 

thac is directly in front of the lens so that its image is 

focused on the center of the plate. The line lying on the 

points P and P„ is the optic axis of the camera. Unless the 
v g 

optic axis is parallel to the z-a::is, the camera will be 
rotated about the optic axis until the image of a vector 
pointing in the +z direction points straight up on the plate. 
If the camera is pointed straight up (or down) , then the 
image of a vector pointing in +x direction will point to the 
right and the image of a vector pointing in the +y direction 
will point up (or down) . 

The plane perpendicular to the optic axis at the lens 
is the principal plane and the plane in which the plate lies 
is the focal plane. The distance between the focal plane and 
the principal plane is the focal length. The focal length 
determines the magnification of the camera. In fact, if a 
segment is perpendicular to the optic axis and if 

H = the length of the segment 

h = the length of the image of the segment 

D = the distance between the segment and the 
principal plane 

d = the focal length of the camera 

then 

h = d*H/D . 
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PTCOMP 


Calling Sequence ; SUBROUTINE PTCOI1P (XI, X2, ISAM) 

Called By ; SHIELD 
II1TCHK 
FINVER 


Local Variables: XI (2): 

Point 1 


X2 (2 ) : 

Point 2 


ISAM: 

(output) ISAM = 0 implies 

XI and 


X2 are different points; 

ISAM = 1 


implies XI and X2 are the 
point. 

same 


Purpose : PTCOMP determines if the points XI and X2 are the 

same point. 
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PTLINE 


Calling Sequence: 

SUBROUTINE PTLINE (XI, X2, XP, IONFLG) 

Called By; I1JTSEC 



Local Variables: 

XI (2) : 

Start of double-precision line 
segment 


X2(2) : 

End of double-precision line 
segment 


XP (2) : 

Double-precision point which may 
or may not lie on the line seg- 
ment 


IONFLG: 

(output) IONFLG = 0 implies XP is 
not on the segment; I01TFLG = 1 
implies XP is on the segment. 


Purpose ; PTLINE determines if point XP (which must lie on 
the line formed by XI and X2) is on the line segment bounded 
by XI and X2. 
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ROUSP 


Calling Sequence : SUBROUTINE ROUSP (IXYZ, IND, IROUS, BUF1, 

IOBPLT, IPER, PER, NX, NY, NZ) 

Called By ; TRILIN 

Local Variables : IXYZ: IXYZ =1, 2 or 3 specifies contours 

through a fixed x, y nr 2 plane 

IND: Fixed value of x, y or 2 ; index in 

grid at which cut is to be taken 

Purpose : ROUSP plots contours of the charge array, ROUS, at 

a specified cut. ROUS is only defined for the inner grid, 
and the 3-D plane at » ich the ROUS values are plotted is 
specified by a fixed x, y or z. IXYZ determines which dimen- 
sion is fixed, and IIID is the fixed value. 

ROUSP also plots the satellite silhouette by reading 
in the building-block data from file IPLT = IABS (IOBPLT) and 
plotting the block perimeters in IPER = PER projected onto 
the appropriate plane. 
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SATPLT 


Calling Sequence : SUBROUTINE SATPLT (NX, NY, NZ, IOBPLT) 

Called By : NASCAP 

Purpose : SATPLT is called when IOBPLT is greater than zero 

indicating satellite illustration plots are desired. 

SATPLT first calls 11ATCAL to set up for and call 14ATPLT 
to generate the material composition satellite silhouettes. 

Next, HIDCEL is called to generate 3-D plots for each 
of three default views, and for any additional views re- 
quested by user input. 


\> 
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SETFAC 


Calling Sequence ; SUBROUTINE SETFAC (I, IDIR, INDFAC, IWORD) 
Called By : SETFWG, SETFOC, SETFTH, SETFFL 

Local Variables : Is I = 1, 2 or 3 implies an x, y, or 

z component of the block face 
normal. 

IDIR: IDIR = 1, 0/ or -1 implies set 

01, 00, or 11 bit pattern for this 
. component. 

INDFAC: INDFAC specifies what face of the 

current block we are defining a 
normal for. 

IWORD: Word currently being filled in with 

normal bit set patterns; following 
completion of this block's face 
normal definitions, IFAC(NBLK) will 
be set equal to IWORD. 

Purpose : The IFAC(80) array stores block face normal informa- 

tion for up to 64 building blocks plus extra face storage for 
up to 16 octagons (which would already be p- t of the set of 
64). For rectangular blocks, the faces are indexed 1 through 
6 and face I will always have the same orientation and normal. 

For blocks which are wedges, octagons, tetrahedra, or 
"FIL111" types, the input block is not symmetric in three 
dimensions, and a given face I of some block may have a normal 
different from the face I normal for a previous block of the 
same type if the two blocks are input with different orienta- 
tions. (The array IPER=PER contains for each block a list of 
vertices, some of which are duplicated, which when traced out 
sequentially in 3-D draws the entire perimeter of the block. 
Face I for a given block type is defined as a particular subset 
of the vertices stored in PER for that block; e.g., vertices 4, 
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3, 8 and 5 night define face 2 of block type 2, vhich is a 
wedge. These face-def inition data are stored in array IV 
in subroutine A2GEN . ) 

A block face normal may be defined by six bits; 0 and 
1, 2 and 3, and 4 and 5 define three bit pairs for the x, y 
and z components of the normal, each pair of which can repre- 
sent either a 1(01), 0(00) or -1(11). All block types except 
the octagon have at most six faces. Consequently, the first 
64 words of IFAC contain, one word (36 bits) per “block, normal 
information for up to six faces per block. The last 16 words 
contain, for a satellite with up to 16 octagons, the normal 
information for the last four octagon faces in bits 0-23. 

For rectangles, because due to their 3-D symmetry of 
general shape, their orientation in the code is always the 
same, the specification of IFAC is a constant set of 36 bits. 
This constant is called IFACR and is defined in RECTAN in a 
data statement. For each rectangular block encountered, 

RECTAN is called and sets IFAC for this block equal to IFACR. 

For the other four block types, however, a somewhat 
complex decoding of input orientation information must be made 
in order to correctly specify the face normal for each block 
in IFAC. For wedges, octagons, tetrahedra and FIL111 type 
blocks, the routines which fill in the set of PER perimeter 
vertices are WEDGES , NICOCT, TETRAH, and FIL111. Each of 
these routines calls, just after PER has been defined, one 
of SETFWG, SETFOC , SETFTH, or SETFFL. Each of these routines 
defines the cell face normal for its block type by filling in 
IFAC(NBLK) with a series of calls to SETFAC. SETFAC fills in 
a particu) ar field of two bits with an 01, 00, or 11, ac- 
cording to whether IDIR equals 1, 0, or -1. Which pair of 
bits of the set of six for a given face is determined by I; 
1=1, 2 or 3 implies the x, y, or z component is being de- 
fined. Which face is being treated (which field of 6 bits) 
is defined by INDFAC. 
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Tx*e building block cell face normals and vertex sets 
are necessary for the definition and sorting of A2 surface 
polygons in subroutine A2GEN. 


SETFFL 


9 


t Calling Sequence ; SUBROUTINE SETFFL 

Called Bv : FIL111 

Purpose ; SETFFL defines the block face normals in array IFAC 
for input blocks which are FIL111 type blocks. A series of 
calls to SETFAC fills in IWDRD, and IFAC (NBLK) is set equal 
to IWORD (See SETFAC). 
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SETFOC 


Calling Sequence : SUBROUTINE SETFOC (N3LK, NOCT, IFAC, II, 

12, 13) 


Called By : NIOOCT 

Local Variables : NBLK: Index of current block being defined. 

NOCT: Number of octagons including current 

one, which have been encountered 
(defines which word of last .16 of 
IFAC will be used in storing normals 
for last 4 faces of this octagon) . 


II: 

12 : 




Some combination of 1, 2, and 3, 
specifying x, y and z components 
of block face normals. 


Purpose : SETFOC defines the block face normals in array IFAC 

for input blocks which are octagons. A series of calls to 
SETFAC fills in IWORD and IFAC (N3LK) is set equal to IWORD 
for the first six faces, while IFAC (64 + NOCT) is set equal 
to a second IWORD for the last four faces (see SETFAC) . 



r. 


SETFTH 

Calling Sequence : SUBROUTINE SETFTH (NBLK, IFAC, ID) 

Called By : TETRAH 

Local Variables ; NBLK: Index of current block being defined. 

ID(3): 1, 0 or -1 for normal specification. 

Purpose : SETFTH defines the block face normals in array IFAC 

for input blocks which are tetrahedra. A series of calls to 
SETFAC fills in IWORD, and IFAC (NBLK) is set equal to IWORD 
(see SETFAC) . 
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SETFWG 


Calling Sequence : SUBROUTINE SETFWG (NBLK , IFAC, ID, II, 12, 

13) 

Called By : WEDGES 

Local Variables: NBLK: Index of current block being defined 


ID (3) : 

1, 0, or -1 for normal 

specification 

11 S ) 

Some combination of 1, 

2, and 3, 

12: [ 

specifying x, y, and z 

components of 

13: ) 

block face normals. 



Purpose : SETFWG defines the block face normals in array 

IFAC for input blocks which are wedges. A series of calls 
to SETFAC fills in IWORD, and IFAC (NBLK) is set equal to 
IWORD (see SETFAC). 
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SHIELD 

Calling Sequence : SUBROUTIIIE SHIELD (HI , N2, Al, A2, ISHLD, 

N3, A3, IMASK, IFF1, IFF2, DV1, DV2, DIR, 
IA1, IA2 , AIR, A2R, IP) 

Called By ; HIDCEL 

Local Variables : Nl: Number of vertices in Al 

N2: Number of vertices in A2 

N3(12): Number of vertices in up to 12 new 

polygons 

Purpose ; SHIELD finishes deciding if an Al polygon is behind 
the current A2 polygon. If so, it masks the Al polygon against 
the A2 polygon, placing new polygons created (representing the 
unshielded portion (s) of the original Al) into the A3 array. 
ISHLD indicates the number of new polygons created by this 
masking. ISHLD = 0 implies Al was completely shielded by A2. 
ISHLD = -1 implies no shielding took place, either because 
there was no overlap or because Al was in front of A2. 

Internal Organization : SHIELD first generates AID and A2D, 

double precision copies of the Al and A2 vertices for use in 
IKSIDD and INTCHK calculations demanding double precision 
accuracy. 

A loop over -the Al vertices determines whether each 
Al vertex is inside A2 or not by calling TNSIDD (on the peri- 
meter is inside) . Within this loop we also perform a behind- 
ness-checking algorithm if IMASK = 0 and if the current vertex 
is inside A2. PRSPPK finds a 3-D point on the line between 
the 3-D surface cell vertex and its 2-D focal plane projection 
(see PRSPPJ) . The 3-D Dine defined by the 3-D lens point (see 
PRSPPJ) and this new 3-D point is the line between the 3-D 
vertex and the 2-D projection. LINPLN recaptures and places 
in XINTA1 the original 3-D coordinates of the Al vertex by 
finding the intersection of the line with the plane of the Al 
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polygon as determined by its normal in IFFI(IAI), bits 0-5, 
and point AIR. AIR is one 3-D vertex of Al; we cannot store 
all 3-D vertices for the Al's due to storage considerations. 
(Note: We do store all 3-D coordinates for the A2 polygons 

in A2R, since there are far fewer of these.) LINPLN is also 
called to find the 3-D point on the plane of A2 which this 
line intersects. This point is placed in XINTA2. Since the 
line from the 3-D vertex through the lens point to the 2-D 
projected point represents the direction of viewing for this 
point, the distances to the viewer of XINTA1 and XINTA2 indi- 
cate whether Al is behind or in front of A2. PRSPPJ calcu- 
lates these two distances, DAI and DA2. 

If the two distances are equal, then this Al vertex is 
actually on A2 in 3-D. We keep track of the number of such 
vertices by incrementing NDA2 for each of these cases. If 
NDA2 reaches 3, indicating that at least throe vertices of 
Al are on A2, the two planes are the same and Al must be on 
A2. We skip out of shield and skip this polygon combination. 
(1'f we masked Al polygons against the A2*s on which they lie, 
the whole satellite would disappear.) 

If DA2 is greater than DAi, then Al is closer to the 
viewer than A2 at this point, which implies that all of Al 
is closer than A2 due to the Al’s being single-cell polygons. 
ISHLD is set to -1 and we skip this combination. 

If DAI is greater than DA2 , then A2 is closer to the 
viewer than Al, and we must mask Al by A2. IMASK is set 
to 1 so ro more behindness testing will be done. 

The inside testing is done for each Al vertex, and tne 
behindness algorithm checked for each Al vertex which is in- 
side A2 untJl IMASK is set to 1. It may happen either that 
no Al vertices are in A2 or that those that are yield only 
inconclusive behindness tests. In this case IMASK may still 
equal zero at the end of this loop. If NDA2 equals 2 at this 
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point# we have Al adjacent to A2 and set ISHLD = -1 to specify 
skipping this combination. If all vertices of Al are inside 
A2, then Al is completely inside A2 (since A2 polygons must 
be convex in this code), and we return with ISHLD = 0. 

A similar loop on A2 vertices follows. For each A2 
vertex# INSIDD determines if the vertex is inside Al. If 
U1ASK = 1# no behindness testing is done. Otherwise, for A2 
vertices inside Al we call LI1JPL2L to find the point on the Al 
plane which the line from the current 3-D A2 vertex to the 
lens point intersects (we need not use PRSPPK to find the 3-D 
vertex since we save all of these for A2's in A2R) . We call 
PRSPPJ to calculate the distance to the viewer of tnis Al 
plane intersection point. We already have this distance for 
the A2 vertex stored in DV2 (another example of information 
stored for A2 polygons but not for Al's). 

We compare these two distances# DA2 and DV2 . If 
DAI = DV2 we increment NDA1. Unless NDA1 reaches three, 
indicating an error, this information is inconclusive.. How- 
ever# if DAI is less than DV2, we may return with ISHLD = -1 
again, since Al is in front of A2; and, if DAI is greater than 
DV2, we set IMASK = 1 to avoid any further behindness tests. 

We determine whether each A2 vertex is inside Al, and 
if IMASK = 0 at this point, check for behindness for A2 
vertices which are in Al# until IMASK is set to 1 or until 
we finish the set of A2 vertices. 

It is possible for IMASK to equal 0 at this point but 
for masking still to be necessary; for example, if two long 
rectangles at right angles cross at their middles, neither 
will have any vertices within the other, hence no behindness 
testing will have been done in the previous two loops. This 
and other somewhat pathological cases must be tested for 
behindness when some intersection point of the two polygons 
is found by the masking algorithm. In this case, an 
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algorithm idertical to that in the first loop is used to 
determine behindness. 

The remainder of the routine masks the A1 polygon by 
the A2 polygon hv following the A1 perimeter around from 
start to finish, making note of all intersections with the 
A2 perimeter, and open: .g and closing new polygons as is 
appropriate. A "new polygon" is a portion of the Al polygon 
which lies completely outside A2 and is composed of one or 
more Al vertices and two intersection points with the A2 
perimeter (these may, of course, also be Al vertices in 
certain cases) . The A2 perimeter intersection points mark 
the start and finish of the "new polygon" at this point, but 
may or may not lie on the same A2 side. If they do lie on 
the same A2 side, then the polygon is already complete. If 
they do not, however, then to complete the new polygon we 
must find that set of A2 vertices which lie on that part of 
the A2 perimeter between the start and finish points of the 
new polygon. Whether the start and finish points lie on the 
same side or not, the task of finishing the vertex set of a 
new polygon is performed by subroutine FINVEk. 

Following the inside tests and behindness checking, a 
loop over the segments of Al begins. A segment of either 
polygon is referenced by the vertex starting it: segment 1 

goes from vertex 1 to 2, segment 2 from vertex 2 to 3, and the 
last segment from vertex IJ to vertex 1. The first vertex in- 
formation is appended to all appropriate lists to facilitate 
handling of the last segment. 

If the last vertex of the current Al segment is out- 
side A2, it is added to the current new polygon vertex list. 
Next, we call INTCHK which checks the current Al segment for 
intersections with all A2 segments. INTCHK returns INT = 
number of distinct intersections with the A2 perimeter. If 
INT =* 0, we go onto the next Al segment. 
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If INT ^ 0, then INTCHK has also returned all inter- 
section points in array XINT, and information (NA2SEG) telling 
how many A2 segments this A1 segment intersected at each 
intersection point {if it intersects at an A2 vertex, then it 
has intersected two A2 segnents at that intersection point). 

It also returns IA2SEG indicating which A2 segment the A1 
segment intersects. FI1JVER will use this information. 

A loop over the number of intersection points follows. 

If the current intersection point is the same as the last 
one treated, we skip it. (There are certain tests here to 
avoid skipping' such a point if it is needed both to start 
the first new polygon and end the last one.) If we keep this 
intersection point, we add it into the current new polygon 
by incrementing INEW (number of points in the new polygon) 
and by moving the intersection point into A3 (INE , NEW) . 

Before we actually move the point into A3, however, we 
must check if, for INEW « 2, we have a segment which is com- 
pletely within A2. If so, we must discard the INEW = 1 point 
and start our new polygon with the current intersection point. 
We determine this by finding the midpoint of the segment 
from the INEW = 1 point tc the INEW = 2 point. Since A2 must 
be convex, if this midpoint is inside A2, we discard this 
segment, while if it is outside, we keep the segment. 

It is at this point that we check for behindness once 
again if IMASK is still equal to zero. 

Since the current point is an intersection point, it 
must either open or close a new polygon. If INEW is equal to 
1, we are opening a new polygon and go on to the next inter- 
section point (if there are no more, we go on to the next 
A1 segment) . 

If INEW is greater than 1, we are closing a new polygon. 
If this is the second intersection point in this new polygon, 
we may call FINVER at this point to finish its vertex set. 
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(The first and last new polygons may need to be added together 
to form a complete new polygon with two intersection points, 
depending on which A1 vertex we started at. ) 

We close this new polygon by incrementing NEW, the 
number of new polygons and by resetting INEW to zero. 

We then check to see if the current intersection point 
is an A2 vertex. If so, we must count it not only as the 
last point of the last new polygon, but also as the opening 
point of the just-opened new polygon. If this is so, we set 
INEW =1, and move the appropriate information into place 
for the new polygon. 

We then loop back to finish any more intersection points 
for the current A1 segment. When these r.re exhausted, we loop 
back to the next A1 segment. When all A1 segment vertices 
and intersection points have been exhausted, we reach state- 
ment label 70. Here we know that there has been at least one 
intersection, and our task ij to determine if either or both 
of the first and last new polygons are incomplete. If the 
first one is incomplete, we join it to the last one. If the 
first one is complete, we either keep or discard the last one 
depending on whether it is complete or not. A polygon to be 
discarded might be, for example, a one -paint polygon. 

If the first and last polygons are added together, 

FINVER is called to complete the vertex set. We are now 
finished and jump to 300 to perform final checks and return 
to HIDCEL. 

If we reach the end of the 200 loop, there have been no 
intersections at all. If this is so, and neither polygon had 
any vertices inside the other, then there is nc overlap and we 
return with ISHLD = -1. 

If all vertices of Al are inside A2, A1 is completely 
shielded and we return with ISHLD = 0 (this has already been 
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checked earlier in the routine) . If all vertices of A2 
are inside A1 and there have been no intersections we error 
return with a RETURN 0. We do_ not allow an A2 polygon to be 
completely enclosed within an Al_ polygon with no intersections . 
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NASCAP COMMON BLOCKS 


Common Block 

Map Segment 

Contents/Purpose 

AN 

FLUXES 

IVN(3,3) ,A (3,3, 3) /Transformation 
matrix for velocities. 

AREAS 

CHARGE 

AREAS (1024) /Sun exposed areas 
of surface cells. 

BLOCKS 

NASCAP 

11(17,33), IBM (33) /Bit arrays. 

CBUF 

TRILIN 

BUF1 (17,17,33) , BUF2 (17 , 17 , 33 ) , 
BUF3 (17 ,17 ,33) /Storage area for 
TRILIN, CHARGE POTENT. 

CD 

MATPLT 

XD, YD. 

CENT 

CHARGE 

ICX (3 ) , CENT (3) /Mesh center co- 
ordinates . 

CIPRT 

POTENT 

IPRT. 

CIT 

NASCAP 

NCYC, MAXITR, INTITR, POTEPS , 
MCYC, IALG. 

CITER 

TRILIN 

ITER. 

CONOPT 

NASCAP 

NCON, ICNOPT (2,8). 

COPT 

NASCAP 

ICREST , IPRi’ST , IOBCAL, ICON, 
ICNVP, I OUTER, ITYPE. 

FLUX 

TRILIN 

Flux definition parameters. 

HIDOPT 

NASCAP 

NDIR, DIKOPT (3,5). 

IGFBLK 

NASCAP 

Graphic output use. 

INDATC 

TRILIN 

ICYC1 , ICYC2 , IG, ISSAVE, NXYZ , 
PCOND f 7) , QCOND (7 ) , QSUMO, IPSAVE 
ICYCS . 

IOCLS 

TRILIN 

IOMAX, IOCLS (128) /Cells for 
output by IOFLX. 

IUNIT 

NASCAP 

Logical unit numbers for external 
files . 

10031 

CNVPLT 

Printer-plotter scratch space. 


Common Block 

Map Segment. 

Contents/Purpose 

LENS 

HIDCEL 

VX, VY, VZ, FL, DFAC. 

LSTSRF 

OBJCOM 

LIST (1024 ) , WPCOND(1024) , W(27,1024: 
Large storage area for OBJDEF. 

MAGDIP 

TRILIN 

BFLDO ( 3 ) , NDIP, PRDIP(6,10)/ 
Magnetic field information. 

MATLS 

NASCAP 

NMAT , MCODE (15 ) , MATPR (20 , 15) / 
Material names and properties. 

MATLS2 

MATPLT 

(identical to MATES) . 

MATOPT 

NASCAP 

NDIV (6) , NDVAL (2,5,6) . 

MCOND 

TRILIN 

CIJ (7,7) , CIJSUM (7 ) , MBIAS, 
VBIAS (7) , MFIX (7) , VFIX(7). 

MESH 

NASCAP 

XMESH, DELTA, ICYC, TIME, DELFAC. 

MISC 

NASCAP 

SWPCND (7) /Potential coefficients. 

MOROPT 

TRILIN 

LONGTM , SHEATH, POTSCL/Keyword 
options . 

NPT 

NASCAP 

NSPTS , NPTS/Number of surface 
and special points. 

NXYZG 

CHARGE 

NX ,NY ,NZ ,NG . 

OUTPOT 

JNOPT 

XS, Y S, ZS. 

PEROBJ 

OBJDEF 

NMCP, PMCND (7, 100) /Multi- 
conductor point potential weights. 

PHCUR 

CHARGE 

IR, 1X1, 1X2, 1X3, JCURCL/Scratch 
variables in photosheath calcu- 
lation. 

PLOT 

CNVPLT 

Scratch area for printer-plotter 
routines. 

PLTOPT 

NASCAP 

ITPART, ITCUR, IROUSP/Plot flags. 

SCRACH 

NASCAP 

Various. 

SRFPRM 

CHARGE 

KTP , ICND, MATNO, MKL ( 3 ) , NCORN , 

' POINTS (3,4) , PCELL, FIELD/Surf ace 
cell parameters. 
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Common Block 


Contents/Purpose 

STOR 

CHARGE 

NBLK, NPBLX (64 ) , P.’QJ (3 , 1024 ) , 
FLOW (1024 ) , FLUXES (1024) , 
SPACER (1024) /Storage area for 
CHARGE segment. 

SUNVEC 

NASCAP 

SUN (3), SFAC/Sun-direction and 
intensity. 

SURF 

OBJ COM 

NSURF , JSURF (1024) /Surface cell 
list. 

SURF 2 

MA'l'PLT 

(Identical to SURF) . 

SURF 3 

CHARGE 

(Identical to SURF) . 

TANK 

TRILIN 

Flux information for test tank 
and reverse trajectory modes. 

TNKSZ 

TRILIN 

IL, ^R, ILGRD/Outer grid 
truncation. 
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11. VARIABLE GLOSSARY 


The following is a glossary of major variables ap- 
pearing in the major subroutines of NASCAP. Subroutines or 
subroutine groups in which a variable is defined or used are 
indicated. Routines in which a variable is merely passed 
along to another level are not normally included. If a 
variable appears in a common block, this is indicated first. 

Variable definitions are as concise and complete as 
possible. However, some reference to documentation of rele- 
vant subroutines or subroutine groups in other sections of 
this report is made where appropriate. 

If a variable .is an input to the code, this is indi- 
cated. Variables in the runstream input file are referred 
to as "user-specified". If a variable is defined in another 
input file such as the keyword input file or the object 
definition file, this is indicated. 
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ALPHA • 


Appears In: CHARGE, POTENT 

Meaning: (CHARGE) Used for potential scaling in 

connection with LONGTIMESTEP option. 

(POTENT) Conjugate gradient algorithm scalar 
equal to r n, r n /u n * (Au) n . 

AREA 

Appears In: CHARGE, HIDCEL 

2 

Meaning: (CHARGE) Area of a surface cell in m . 

AREAS; AREAS (1024) 

Appears In: CHARGE 

Meaning: Common block containing array of same name 

containing fraction of each surface cell 
exposed to sunlight. 

AREA1 (1024) 

Appears In: HIDCEL 

Meaning: First used to hold final shielded 2-D areas 

of surface cells; then holds those areas 
divided by original unshielded 2-D areas 
to yield final fractional areas with the 
dot product (angle) effect removed. 

AREA2 (1024) 

Appears In: HIDCEL 

Meaning: Original unshielded 2-D areas of surface 

cells. 

AUNCON ( 7 ) 

Appears In: POTENT and related routines. 

Meaning: Conjugate gradient cc efficient product (Au) n 

for conductors. 
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Al (2,12 , 600) 


Appears In: HIDCEL and related routines. 

Meaning: A1(1,I,J) and A1(2,I,J) are x and y co- 

ordinates of Ith vertex of Jth A1 polygon. 


Appears In: HIDCEL and related routines. 

Meaning: A1B(I) is distance to viewer of vertex 

furthest from viewer for the surface cell 
which projects into the Jth A1 polygon. 

Appears In: HIDCEL, ADDA1 

Meaning: Temporary save of AlB(IAl). 


Appears In: HIDCEL 

Meaning: Temporary save of AlB(IAl). 

!L 

Appears In; HIDCEL and related routines. 

Meaning: A1F(I) is distance to viewer of vertex 

closest to viewer for the surface cell which 
projects into the Ith A1 polygon. 


A1FT 


Appears In: HIDCEL, ADDA1 

Meaning: Temporary save of AlF(IAl). 


Appears In: HIDCEL 

Meaning: Temporary save of AlF(IAl). 
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A1R(3, 600) 


Appears 
Meaning : 


AIRS (3) 

Appears 
Meaning : 

A2 (2,9,150) 

Appears 
Meaning : 

A2B (150 ) 

Appears 

Meaning: 


A2BT 

Appears 
Meaning : 

A2F (150) 

Appears 
Meaning : 


A2FT 

Appears 
Meaning : 
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In: HIDCEL and related routines 

A1R(1,I), AIR (2 , 1 ) , and A1R(3,I) are x, y 
and z coordinates of first vertex of 3-D 
surface cell which projects into Ith Al 
polygon. Used with surface cell normal to 
define plane of cell for LINPLN (see 
SHIELD) . 

In: HIDCEL, ADDA1 

Temporary save of (AIR (K, IA1) ,K=1, 3) . 

In: HIDCEL and related routines. 

A2 (1, 1, J) and A2(2,I,J) are x and y co- 
ordinates of Ith vertex of Jth A2 polygon. 


In: HIDCEL and related routines. 

A2B (I) is distance to viewer of vertex 
furthest from viewer for the building-block 
face which projects into the Ith A2 polygon 


In: HIDCEL 

Temporary save of A23(IA2). 


In: HIDCEL and related routines. 

A2F(I) is distance to viewer of vertex 
closest to viewer for the building-block 
face which projects into the Ith A2 polygon 

In: HIDCEL 

Temporary save of A2F(IA2). 
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A2R (3,9,150) 


Appears In: HIDCEL and related routines. 

Meaning: A2R (1 , 1, J) ,- A2R (2, 1 , J) , A2R(3,I,J) are x, y 
and z coordinates of Ith vertex of building- 
block face which projects into Jth A2 polygon. 

A3 (2,24,12) 

Appears In: HIDCEL and related routines. 

Meaning: A3(1,I,J), A3(2,I,J) are x and " coordinates 

of Ith vertex of Jth new polygon being 
created by SHIELD. 

BETA 

Appears In: CHARGE, POTENT 

Meaning: (CHARGE) Factor used in LONGTIMESTEP 

correction. 

(POTENT) Conjugate gradient algorithm scalar 
equal to -r n * (Au) n /u n * (Au) n . 

BUF1 (17,17,33) 

Appears In: C0MM0N/C3UF/; POTENT and related routines; 

other routines. 

Meaning: First of three identical buffer storage arrays 

used mainly for conjugate gradient array 
storage of one grid's worth of data at a tine. 
Used in other parts of code for other pur- 
poses as well, sometimes under different names. 
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BUF2 (17,17,33) 


Appears In: COMMON/CBUF/; POTENT and related routines; 

other routines. 

Meaning: Second of three identical buffer storage ar- 

rays used mainly for conjugate gradient array 
storage of one grid's worth of data at a time. 
Used in other parts of code for other pur- 
poses as well, sometimes under different names. 

BPF3 (17,17,33' . 

Appears In: COMMON/CBUF/; POTENT and related routines; 

other routines. 

Meaning: Third of three identical buffer storage arrays 

used mainly for conjugate gradient array 
storage of one grid's worth of data at a time. 
Used in other parts of the code for other pur- 
poses, sometimes under different names. 

C 

Appears In: POTENT and related routines. 

Meaning: Conjugate gradient algorithm scalar (for 

algorithm 1): C 1+1 = 1. + BETA*C 1 . 

CENT; C E NT (3) 

Appears In: CHARGE 

Meaning: Common block containing array with values 

(NX/2+1), (NY/2+1), (NZ/2+1 ) . 

CIJ (7,7) 

Appears In: CCMMON/MCOND/ ; POTENT and related rou- 

tines; INKEYW. 

Meaning: CIJ(J,J) is capacitance between Ith and Jth 

conductors associated with glue joint or 
insulating spacer. 
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CIJSUM (7 ) 


CNET 


CNET2 


CS 


DELFAC 


DELTA 


DFAC 


Appears In: COMMON/MCOND/; POTENT and related routines 

INKEYW. 

Meaning: CIJSUM(I) = 2 CIJ(I,J). 

J 


Appears In: CHARGE 

Meaning: Net charging current. 


Appears In : CHARGE 

Meaning: Net charging current at advance timestep 

computed by LONGTIMESTEP option. 


Appears In: POTENT 

Meaning: Temporary save of C 1 after C has been set 

1 i i 


Appears In: CHARGE 

Meaning: Factor by which DELTA is multiplied after 

each cycle; runstream input. 

Appears In: CHARGE 

Meaning: Length of timestep in seconds; runstream 

input. 


Appears In: HIDCEL 

Meaning: Constant specifying distance to viewer from 

satellite center. 
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DIR(3,3) 


Appears In: SATPLT, HIDCEL 

Meaning: (SATPLT) DIR(1,I), DIR(2,I) andDIR(3,I) 

are x, y and z components of the direction- 
of-view vector pointing from the center of 
the satellite toward the viewer, for the 
Ith default satellite illustration view. 

(HIDCEL) DIR (1) , DIR (2) , and DIR{3) are 
the x, y and z components of the normalized 
direction-of-view vector for the current 
HIDCEL cell. 

DIROPT (3 , 5 ) 

Appears In: SATPLT 

Meaning: DIR(1,I) , DIR(2,I), and riR(3,I) are the 

x, y and z components of the Ith user- 
specified optional satellite illustration 
view. 

DNORM 

Appears 'In: HIDCEL 

Meaning: Length of unnormal j zed direction-of-view 

vector PDIR (3) input to HIDCEL. 

DUM 

Appears In: NASCAP, 3HARGE, HIDCEL, CONPLT, ROUSP, 

PARPLT 

Meaning: Dummy read variable for quick skip of a 

binary record on an external file. 

DV1 

Appears In: HIDCEL, A1GEN 

Meaning: Not currently used. 
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PV2 (9,150) 

Appears In: HIDCEL and related routines. 

Meaning: DV2(I,J) is the distance to viewer of the 

Ith vertex of the building-block face 
which projects into the Jth A2 polygon. 


FLOW (1024) 

Appears In: CHARGE 

Meaning: Low energy electron flux emitted from each 

surface cell. 


FLUKES (1024) 

Appears In : CHARGE 

Meaning: Net flux to each surface cell. 


IALG 

Appears In: INOPT, TRILIN, POTENT 

Meaning: User-specified input determining which 

conjugate gradient potential solution 
algorithm will be used: IALG = 0 im- 

plies algorithm 0; IALG = 1 implies 
algorithm 1. 
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IAREA 


Appears In: COMMON/IUNIT/,* NASCAP, INOPT, HIDCEL 

Meaning: User-specified I/O unit number of external 

file on which shadowed surface cell areas 
(AREAl) are written. If IAREA > 0, a new 
sun direction vector (SUN) is used to make 
a new shadowing calculation. If IAREA < 0, 
the old shadowed areas on unit IABS (IAREA) 
are used. When HIDCEL is called from 
SA'fPLT with argument IAREA = 0 , no area 
calculations are made, and the non-hidden 
line building -block plot is made in addi- 
tion to the hidden line cell-by-cell plot. 

IAUN 

Appears In: COMMON/ INN IT/; INOPT; POTENT and related 

routines. 

Meaning: User-specified I/O unit number of external 

file containing (Au) n , conjugate gradient 
coefficient product values for grid points. 

IA1 

Appears In: HIDCEL and related routines. 

Meaning: Index of current A1 polygon. 

IA1T 

Appears In: HIDCEL 

Meaning: Index in original surface cell list of 1024 

cells of the surface cell which generated 
parent (original, ur. shielded) A1 polygon of 
current A1 polygon, where final polygon in- 
formation is printed. 
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IA1TT 


Appears In: HIOCEL 

Meaning: Index in original surface cell list of 1024 

cells of the surface cell which generated 
parent (original, unshielded) A1 polygon of 
current A1 polygon in shielding loops. 

IA2 

Appears In: HIDCEL and related routines. 

Meaning: Index of current A2 polygon. 

IBM (33) 

Appears In: NASCAP, CHARGE, OB JDEF , /BLOCKS/ 

Meaning: A bit mask array: IBM(I) = 2** (1-1) . 

IC 

Appears In: HIDCEL, OB JDEF, POTENT 

Meaning: (HIDCEL) Loop index over surface cell area 

list; 1-1024. 

(OB JDEF ) Loop index over conductors; 1-7. 
(POTENT) Loop index over conductors; 1-7. 

ICNV P 

Appears In: COMMON/COPT/ ; INOPT, TRILIN 

Meaning: User-specified input flag indicating 

whether or not and how often potential 
convergence printer plots are desired: 

ICNVP <_ 0 implies no plots; ICNVP > 0 
implies plot convergence data every 
ICNVP time cycles, starting with 
cycle 1. 


ICON 


Appears In: COMMON/COPT/; NASCAP, INOPT, TRILIN 

Meaning: User-specified input flag indicating whether 

or not and how often potential contour plots 
are desired. ICON <_ 0 implies no plots; 

ICON > C implies plot potential contours 
every ICON time cycles, starting with 
cycle 1. 

ICREST 

Appears In: COMMON/COPT/; INOPT, INDATA, TRILIN 

Meaning: User-specified new time cycle restart flag: 

ICREST = 0 implies fresh start of new 
problem if IPREST = 0, or potential cal- 
culation restart at same time cycle as 
last run if IPREST / 0. ICREST = 1 im- 
plies restart at beginning of the time 
cycle following last time cycle of pre- 
vious run. 

ICYC 

Appears In: COMMON/MESH/; TRILIN, CHARGE, POTENT 

Meaning: Index of current time cycle. 

ICYCS 

Appears In: COMMON/ INDATC/; INDATA, TRILIN 

Meaning: Temporary save of last cycle completed. 

Saved upon run completion for restarts. 

ICYC1 

Appears In: COMMON/ INDATC/; INDATA, TRILIN 

Meaning: Index of time cycle at which this run will 

begin. Equals 1 for a fresh start. Equals 
(ICYC2 of last run + 1) ' for' T’c^f# ij re start . 
Equals ICYC2 of last run for a potential 
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ICYC2 


Appears In: COMMON/ INDATC/ ; INDATA, TRILIN 

Meaning: Index of last time cycle which this run 

will complete. Equals ICYC1 + NCYC - 1. 

I FAC (80) 

Appears In: COMMON/LSTSRF/ (in OBJDEF routines); 

RECTAN, WEDGES, NIOOCT, TETRAH, FIL111, 
SETFWG, SETFOC, SETFTH, SETFFL, HIDCEL, 
A2GEN 

Meaning: IFAC(I) = 1,64 contains bit encoded data 

specifying the block face normals for the 
first six faces of the Ith satellite build- 
ing block. The 36 bits form six groups of 
six bits each; for a group of six, the 
first two bits specify the x component, 
the second two bits the y component, and 
the third two bits the z component of the 
normal as follows: 00 = 0; 01 = 1; 

11 = -1; 10 not permitted. The last 
sixteen words of IFAC contain bit encoded 
data for normals for the last four faces 
of up to sixteen octagonal blocks (the 
only block type with more than six faces) . 
See subroutine SETFAC. 


IFFS 

Appears In: HIDCEL, ADDA1 

Meaning: Temporary save of I7F1(IA1). 
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IFFM600) 

Appea.es In: HIDCEL and related routines. 

Meaning: IFFl(I) contains bit encoded data indi- 

cating the surface cell normal and index 
for the surface cell which projected into 
the Ith II polygon. Bits 0,1: x component 

of normal; Bits 2,3: y component of nor- 
mal; Bits 4,5: z component of normal, 

where 00 = 0, 01=1, 11 = -1, 10 not 
permitted. Bits 6-15: index (1-1024) 

in surface cell list of parent surface 
cell. 

IFF2 (150) 

Appears In: HIDCEL and related routines. 

Meaning: The first six bits of IFF2 (I) contain bit 

encoded data indicating the block face nor- 
mal of the building-block face which pro- 
jected into the Ith A2 polygon, exactly as 
IFF1 stored A1 polygon normals. 

II (.17,33) 

Appears In: NASCAP, OBJDEF, OB JSPC, /BLOCKS/ 

AND (II (IY, IZ) , IBM (IX)) / 0 if 
point (IX , IY , IZ ) is a special, sur- 
face or interior point. 

IM 

Appears In: POTENT 

Meening Loop index over biased conductors in 
ROUSCN definition. 
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IMASK 


Appears In: HIDCEL, SHIELD 

Meaning: IMASK = 0 implies we do not yet know if 

A1 is behind A2. IMASK = 1 implies we 
have ascertained that A1 is behind A2 
and we must therefore mask A1 by A2. 

See subroutines HIDCEL and SHIELD. 

INTITR 

Appears In: COMMON/CIT/; INOPT, POTENT 

Meaning: User-specified input flag specifying 

that the conjugate gradient algorithm 
is to be reinitialized every INTITR 
iterations (not normally used) . 

IOBCAL 

Appears In: COMMON/COPT/; INOPT, NASCAP 

Meaning: User-specified i-put flag for OBJDEF call: 

IOBCAL = 1 implies call OBJDEF and proceed 
normally. IOBCAL = 0 implies skip OBJDEF 
call and proceed normally (in this case 
OBJDEF must have been called for the cur- 
rent satellite in a previous run) . 

IOBCAL = -1 implies call OBJDEF (and 
SATFLT and HIDCEL if called for) , and 
stop. 

IOBJ 

Appears In: COMMON/ IUN IT/; INOPT, OBJDEF, OBJSPC 

Meaning: I/O unit number of external file on which 

OBJDEF writes the PMCND, LIST, WPCOND, and 
W arrays for use by OBJSPC. 
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IOBPLT 


Appears In: COMMON/IUNIT/; INOPT, OBJDEF, NASCAP , 

ETNGUN , SOURCE/ CONPLT, PARPLT , ROUSP, 
HIDCEL/ MATPLT 

Meaning: User-specified I/O unit number of external 

file on which OBJDEF writes satellite de- 
finition data, and from which all plot rou- 
tines, etc., reading this data read it. 

Also a flag indicating whether or not satel- 
lite illustration plots are desired. 

IOBPLT > 0 implies plot material composition 
satellite silhouettes, non-hidden line satel- 
lite building-block plots, and hidden-line 
cell-by-cell plots. IOBPLT < 0 implies do 
not plot, and use unit number IPLOT = 

IABS (IOBPLT) . 

IP 

Appears In: COMMON/IUNIT/; INOPT, TRILIN, POTENT and 

related routines, CHARGE and related rou- 
tines, CONPLT 

Meaning: User-specified I/O unit number of external 

file in which grid potentials are stored. 

IPART 

Appears In: COMMON/IUNIT/; INOPT, NASCAP, TRILIN, PUSH, 

PARPLT 

Meaning: User-specified I/O unit number of external 

file cn which PUSH writes particle trajectory 
data for plotting by PARPLT. Also a flag 
indicating wnether or not particle trajectory 
plots are desired. IPART > 0 implies write 
particle trajectory d=ta for particles landing 
in cell I0CL(1) onto file IPART, and plot 
these trajectories for the first cycle of 

314 the run in subroutine PARPLT. 


IPER(3,32,64) 


Appears In: COMMON/LSTSRF/ (in OBJDEF and related rou- 

tines); HIDCEL, CONPLT , RQUSP, PARPLT, 
ETNGUN, SOURCE 

Meaning: Same storage as PER in all routines except 

OBJDEF. Integer version of data in PER. 
IPER(1,I,J), IPER(2,I,J), and IPER(3,I,J) 
are the x, y and z coordinates of the Ith 
vertex in the complete-perimeter vertex 
list for satellite block J. For each block, 
OBJDEF generates a complete-perimeter vertex 
list, containing some duplicated vertices, 
which when traced out in sequence defines 
the entire block perimeter. 

IPFLAG 

Appears In: NASCAP 

Meaning: IPFLAG is set to one if any plots are called 

for by user input; otherwise it is zero. 
IPFLAG is used to determine whether a plot 
file must be opened and closed. 

IPLOT 

Appears In: NASCAP, INOPT, OBJDEF, CHARGE, ETNGUH, 

SOURCE 

Meaning: IPLOT = IADS (IOBPLT) . 
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IPREST 


Appears In: COMMON/COPT/; INOPT, TRILIN, INDATA 

Meaning: User-specified potential restart flag. 

IPREST = 0 implies begin with new poten- 
tial calculation at next cycle (cycle 1 if 
ICKEST = 0, or cycle following last cycle 
of previous run if ICREST = 1) . 

IPREST = 1 implies restart in the middle 
of the previous run's last cycle potential 
calculation by re-initializing the conjugate 
gradient algorithm with the last run's final 
potentials as the initial guess. Take the 
last run's potentials from file IP. 

IPREST < 0 implies restart as for IPREST = 1, 
but take the last run’s potentials from file 
IPRDP = IABS (IPREST) . 

IR 

Appears In: COMMON/ IUN IT/; INOPT; POTENT and related 

routines. 

Meaning: User-specified I/O unit number of external 

file on which the conjugate gradient array 
r (for grid points) is stored. 

IREMOV (600) 

Appears In: HIDCEL, A1COMP 

Meaning: A1 polygon removal flag. If A1 polygon IA1 

is partially or wholly shielded by the cur- 
rent A2 polygon, then IREMOV(IAl) is set to 1. 
Otherwise IREMOV (IA1) = 0. 

IREMOV is used by A1COMP to remove unneces- 
sary A1 polygons from the A1 information 
arrays. See subroutines HIDCEL and A1COMP. 
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IROUS 


Appears Ins COMMON/IUNIT/; INOPT, POTENT and related 
routines; CHARGE and related routines. 

Meaning: User-specified I/O unit number of external 

file on which the ROUS charge array lies. 

IROUSP 

Appears In: CGMMON/PLTOPT/; INOPT, NASCAP, TRILIN 

Meaning: User-specified plot flag for ROUS contour 

plots. IROUSP =0 implies no plots. 

IROUSP / 0 implies plot final ROUS contours. 

ISAT 

Appears In: COMMON/IUNIT/; INOPT, OBJ DEP and related 

routines . 

Meaning: User-specified I/O unit number of the ex- 

ternal file containing satellite definition 
input to OBJDEF, 

ISHLD 

Appears In: HIDCEL, SHIELD 

Meaning: ISHLD _> 0 implies rhat masking of the cur- 

rent Al polygon by the current A2 polygon 
has produced ISHLD new polygons (if 
ISHLD = 0 the Al polygon was completely 
masked by the A2 polygon). ISHLD = -1 
implies no masking occurred either because 
Al was in front of A2 or because the poly- 
gons did net overlap. 

ISPARE 

Appears In: COMMON/IUNIT/; INOPT, POTENT and related 

routines . 

Meaning: I/O unit number of the external file *,nich 

is currently "spare", i.e., whose data is 
currently obsolete. See description of con- 
jugate gradient potential solution. 
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ITCUR 


Appears In: COMMON/PLTOPT/ ; INOPT, NASCAP, SOURCE 

Meaning: User-specified plot flag for tank test cur- 

rent contours. ITCUR = 0 implies no plots. 
ITCUR / 0 implies plot current contours. 

ITER 

Appears In: COMMON/CITER/; POTENT and related rou- 

tines. 

Meaning: Index of current potential iteration in 

conjugate gradient solution. 

ITPART 

Appears In: COMMON /FLTOPT/ ; INQPT , NASCAP, ETNGUN 

Meaning: User-specified plot flag for tank test 

particle trajectory plots. ITPART = 0 
implies no plots. ITPART / 0 implies 
plot particle trajectories. 

ITYPE 

Appears In: COMMON/COPT/ ; INOPT, CHARGE, FLXDEF 

Meaning: User-specified problem type flag. 

ITYPE = 1 implies test tank case. 

ITYPE = 2 implies Maxwellian isotropic 
plasma approximation. 

ITYPE = 3 implies particle-pushing space 
case with experimental (ATS-5-) 
ambient plasma. 

IU 

Appears In: COMMON/IUNIT/; INOPT, POTENT and re- 

lated routines. 

Meaning: User-specified I/O unit number of external 

file on which the conjugate gradient array 
u (for grid points) is stored. 
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IV 


Appears In: HIDCEL 

Meaning: Loop index over vertices of a final A1 poly- 

gon when finding 2-D minima and maxima for 
plot user area definition. 

IWARN 

Appears In: POTENT, TRILIN 

Meaning: (for S-Cubed only) At the start of the 

potential iteration loop, S3 WARN is called 
to determine if the number used up all but 
a pre-sp*ecified amount of its allotted cpu 
time. If it has, IWARN is seu to one which 
tells the run to attend to necessary plots 
and restart dumps, and then terminate. If 
not, IWARN remains zero and the run proceeds 
normally. 

IY 

Appears In: COMMON/ IUH IT/; INOPT, POTENT and related 

routines. 

Meaning: User-specified I/O unit number of external 

file on which the conjugate gradient array 
y (for grid points) is stored. 

11 

Appears In: HIDCEL 

Meaning: Index of new polygons created by last SHIELD 

call running from 1 to ISHLD as the new poly- 
gons are added to A1 by ADDA1 . 

JC 

Appears In: POTENT 

Meaning: Inner loop index (1 to MBIAS) of double loop 

over biased conductors which defines SWPCIJ 
and SSCIJ . 
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JM 


Appears In: POTENT 

Meaning: Inner loop index (1 to MBIAS) of double 

loop over biased conductors which defines 
ROUSCN. 

JSURF (1024 ) 

Appears In: CHARGE, HIDCEL, OBJDEF 

Meaning: An array containing the internal codes 

for surface cells (Figure 6.14). 


MAXCYC 


Appears In: POTENT 

Meaning: Parameter equal to 100 which defines size 

of RITER, UITER and PCITER storing poten- 
tial convergence printer plot data for up 
to 100 potential iterations. 

MAXIDF 

Appears In: INOPT 

Meaning: Default value of MAXITR. If MAXITR as 

read in is negative, MAXITR will be set 
equal to MAXIDF = 50. 

MAXITR 

Appears In: COMMON/ C IT/ ; INOPT, POTENT 

Meaning: User-specified limit on the number of 

potential iterations. Currently there 
are no convergence tests and POTENT 
will iterate MAXITR times no matter 
what (unless S3WARN indicates cpu time 
is running low) . 
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MB IAS 


MBIAS1 


MCYC 


MFl'X (7) 


M1NITR 
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Appears In: COMMON/MCOND/ ; INKEYW, POTENT and 

related routines. 

Meaning: Index of last biased conductor (= number 

of biased conductors including conductor 
one). Defined via the Keyword Input file. 

Appears In: POTENT 

Meaning: MBIAS + 1. Index of first r.on-biased 

conductor. 

Appears In: COMMON/CIT/ 

Meaning: User-specified flag indicating whether or 

not and how often the potential solution 
will be performed. MCYC = 0 implies no 
potential solution; old or scaled potentials 
will be used. MCYC > 0 implies potentials 
will be solved when ICYC = 1 mod MCYC. 


Appears In: COMMON/MCOND/; INKFYW, POTENT and 

related routines. 

Meaning: Fixed potential conductor flag. 

MFIX(I) = 0 implies conductor I is not 
fixed. MFIX{I) = 1 implies conductor I 
is fixed. 

Appears In: POTENT 

Meaning : Iteration at which the potential iteration 

loop will begin. Normally equals one, but 
if INTITR is non-zero, MINITR may be reset 
to (N* INTITR) + 1. 
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NA1 


NAIM 


NA1HEW 


NA? 


NBLK 


Appears In: HIDCEL and related routines. 

Meaning: Number of A1 polygons currently in the A1 

information arrays. 


Appears Ir.: HIDCEL and related routines. 

Meaning: Parameter equal to 600 used in defining A1 

information array sizes. No more than NAIM 
A1 polygons may exist before list compression 
by A1COMP. 

Appears In: HIDCEL, A1COMP 

Meaning: Output from ADDA1 indicating the total num- 

ber of old A1 polygons plus the number of new 
ones which have been appended to the A1 lists 
by masking with the current A2 polygon. See 
HIDCEL. 


Appears In: HIDCEL and related routines. 

Meaning: Total number of A2 polygons. 


Appears In; C0M1I0N/LSTSRF/ (in OBTDEF) ; 03JDEF , HIDCEL, 

A2GEN, CONPLT, ROUSP, PARPLT, ETNGUN, SOURCE 
Meaning: Number of satellite definition blocks 

("fcyilding-blocks") as defined in the OEJDEF 
input file ISAT. 
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Appears In: OBJDEF, POTENT and related routines; 


NCMAX 


NCON 


NDIR 


NG 


NITERP 


Meaning: Number of conductors in current run. 

Appears In: NASCAP, INOPT, OBJDEF and related routines; 

TRILIN, INDATA, CHARGE and related routines 
POTENT and related routines. 

Meaning: Parameter equal to 7 . Maximum allowable 

number of conductors. 

Appears In: COMMON/C ONOPT/; INOPT, TRILIN 

Meaning: User-specified number of optional potential 

contour plot cuts desired (normally zero) . 

Appears In: COMMON/HIDOPT/; INOPT, SATPLT 

Meaning: User-specified number of optional satellite 

illustration views for the three-dimensional 
plots. 


Appears In: Entire code. COMMON/NXYZG/ (in CHARGE and 

related routines). 

Meaning: User-specified number of nested NX x NY x NZ 

grids in the problem. 

Appears In: TRILIN, POTENT 

Meaning: Number of potential iterations for which 

potential convergence printer plot data has 
been stored. 



NMCMAX 


Appears In: OBJDEF and related routines; POTENT and 

related routines. 

Meaning: Parameter equal to 100. Maximum allowable 

number of multiple-conductor points. See 
OBJDEF. 

NPBLK (64 ) 

Appears In: COMMON/LSTSRF/ (in OBJDEF); OBJDEF and 

related routines; HIDCEL, A2GEN, CONPLT, 
ROUSP, PARPLT, ETNGUN, SOURCE. 

Meaning: NPBLK (I) is the number of points in the 

IPER complete perimeter list for satel- 
lite input block I. See IPER. 

NPTS 

Appears In: NASCAP, OBJDEF 

Meaning: The number of special points. 

NSPTS 

Appears In: NASCAP, OBJDEF 

Meaning: The number of surface points. 

NSURF 

Appears In: CHARGE, HIDCEL, OBJDEF 

Meaning: The number of surface cells. 

NTOTMX 

Appears In: OBJDEF and related routines; OBJSPC. 

Meaning: Parameter equal to 1024. Maximum allowable 

number of surface cells, surface poirts, or 
special points. 

NV 

Appears In: HIDCEL and related routines. 

Meaning: Temporary storage of NV1(I) for use as a 

do loop limit. 
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NV1 (600) 


Appears In: HIDCEL and related routines. 

Meaning: HV1(I) is the number of vertices in the Ith 

A1 polygon (including the repeat of the first 
vertex at the end of the list) . See HIDCEL. 

N V2 (150) 

Appears In: HIDCEL and related routines. 

Meaning: NV2(I) is the number of vertices in the Ith 

A2 polygon (including the repeat of the first 
vertex at the end of the list) . See HIDCEL. 

N VER 

Appears In: HIDCEL and related routines. 

Meaning: Parameter equal to 12. Maximum allowable 

number of vertices in any A1 polygon (in- 
cluding the repeat of the first vertex at 
the end of the list) . Although all A1 
polygons start out with a maximum of five 
vertices, as they are masked, they can be- 
come concave and/or many-sided. If the code 
attempts to generate a polygon with twelve 
distinct vertices (which will need thirteen 
spaces in the vertex array) , it will error 
terminate. A2 polygons are static and have 
a maximum of nine vertices (3 + 1) . 


KVM1 

Appears In: HIDCEL 

Meaning: Number of vertices minus one for some poly- 

gon (hence the number of distinct vertices) . 

NVP1 

Appears In: HIDCEL 

Meaning: Number of vertices plus one for some polygon. 
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NX 


NXM 


NY 


NYM 


NZ 


N3K 


Appears In: Entire code. COMIiON/NXYZG/ (in CHARGE 

and related routines). 

Meaning: User-specified number of x direction grid 

points in each of the nested grids; 
normally 17. 


Appears In: OBJDEF and related routines; TRILIN. 

Meaning: Parameter equal to 17. Maximum allowable 

number of x direction grid points. 


Appears In: Entire code. COMMON/NXYZG/ (in CHARGE 

and related routines) . 

Meaning: User-specified number of y direction grid 

points in each of the nested grids; 
normally 17. 

Appears In: OBJDEF and related routines; TRILIN. 

Meaning: Parameter equal to 17. Maximum allowable 

number of y direction grid points. 


Appears In: Entire code. COMMON/NXYZG/ (in CHARGE 

and related routines). 

Meaning: User-specified number of z direction grid 

points in each of the nested grids; an 
integer of form 4n+l, 4 <_ n <_ 8. 

Appears In: OBJDEF and related routines; TRILIN. 

Meaning: Parameter equal to 33. Maximum allowable 

number of z direction grid points. 
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N3 (12) 


t 
i 

Appears In: HIDCEL, SHIELD and related routines. 

Meaning: N3(I) is the number of vertices in the 

Ith new polygon being created by SHIELD. 

See subroutines HIDCEL and SHIELD. 

N3T 

Appears In : HIDCEL 

Meaning: Temporary save of N3 (I) for use as a do 

loop limit. 

PC ITER (100) 

Appears In: POTENT, TRILIN 

Meaning: Potential convergence printer plot storage 

for up to 100 iterations* worth of the 
potential on conductor one. 

PCOND ( 7 ) 

Appears In: POTENT and related routines; TRILIN, 

INDATA, CHARGE and related routines. 

Meaning: PCOND (I) is the potential on conductor I. 

PDIR(3) 

Appears In: HIDCEL 

Meaning: Unnormalized direction-of-view vector for 

this call to HIDCEL. See DIR(3). 

PEPSDF 

Appears In: INOPT 

Meaning: Designed as a default epsilon for poten- 

tial convergence tests. If POTEPS is read 
in as 0, it is set equal to PEPSDF. At 
present, no potential convergence tests 
are included in the code hence neither 
PEPSDF or POTEPS is used. 



PER (3 ,32, 64) [or PER(32,64) in OBJDEF] 


Appears In: COMMON/ JSURF/ (in OBJDEF and related 

routines); HIDCEL, CONPLT, ROUSP, 
PARPLT, ETNGUN, SOURCE. 

Meaning: In all routines except OBJDEF and its 

descendants, PER is a real copy of the 
integer information in IPER (in fact, 
it is the same storage area). In 
OBJDEF, PER holds the same information, 
but the data is bit-encoded three 
integers to a word: Bits 5 to 10: 

x coordinate; Bits 11 to 16: y co- 
ordinate; Bits 17 to 22: 2 coordinate. 

See IBM (33). 

PHOJ (3,1024 ) 

Appears In: CHARGE 

Meaning: An array containing the photosheath cur- 

rents in volume cells immediately sur- 
rounding the object. 


PI 

Appears In: POTENT 

Meaning: 3.1415926. 

POINTS (3,4) 

Appears In: CHARGE 

Meaning: An integer array containing the vertices 

of current surface cell. 

POTSCL 

Appears In: TRILIN 

Meaning: A variable containing the potential 

scaling option literal . 
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P3 (36 , 36 ) 


Appears In: POTENT and related routines. 

Meaning: Scratch storage for the I-prefix and 

O-prefix routines called by subrov tines 
IRC ALL and ORCALL. 


Appears In: CHARGE and related routines; TRILIN, 

INDATA, POTENT and related routines. 
Meaning: QCOND(I) is the charge on conductor I. 


Appears In: HIDCEL, PRSPPJ 

Meaning: Real 3-D x coordinate of center of 

inner grid. Used by subroutine PRSPPJ 
(entry point PPJSET) . 


Appears In: HIDCEL, PRSPPJ 

Meaning: Real 3-D y coordinate of center of 

inner grid. Used by subroutine PRSPPJ 
(entry point PPJSET). 


Appears In: HIDCEL, PRSPPJ 

Meaning: Real 3-D z coordinate of center of 

inner grid. Used by subroutine PRSPFJ 
(entry point PPJSET) . 


Appears In: TRILIN 

Meaning: The total charge in the inner grid. 


Appears In: TRILIN 

Meaning: QSUM — (sum of charges on fixed potential 

conductors) . 


RCOND ( 7 ) 


Appears In: POTENT and related routines. 

Meaning: Conjugate gradient algorithm r array 

values for each of up to seven con- 
ductors. See subroutine POTENT. 


RDOTAU (double precision) 


Appears In: POTENT, RUPDAT 

Meaning: Conjugate gradient array dot product 

r n+ * # (Au) n used to calculate BETA. 
See subroutine POTENT. 

RDOTR (double precision) 

Appears In: POTENT, URSETO, RUPDAT 

Meaning: Conjugate gradient array dot product 


RDOTRS (double precision) 

Appears In: POTENT, RUPDAT 

Meaning: Temporary save of previous iterations 

RDOTR after RDOTR has been set to 
r n+l. c n+l_ 

RDOTS 

Appears In: POTENT 

Meaning: Temporary save of first iteration RDOTR 

for possible convergence testing (not 
implemented) . 

RITER(IOO) 


Appears In: POTENT, TRILIN 

Meaning: Potential convergence printer plot 

storage for up to 100 iterations' worth 
of RDOTR. 
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ROUS (17/17,33) 


In: 


TRILIN, INDATA, CHARGE and related 
routines; POTENT and related rou- 
tines. 


i. 


Appears 


Meaning: Amounts of charge at inner grid points. 

Right-hand-side of potential equations. 
This data is also stored in file IROUS. 
See subroutine CHARGE. 

ROUSCN (7) 

Appears In: POTENT and related routines. 

Meaning: Right-hand-side of linear equations for 

the conductor potentials. 


RTDOT 

Appears In: POTENT 

Meaning: Monotonic residual dot product with it- 

self, r*r, created by potential algo- 
rithm one for judging quality of poten- 
tial convergence. 


SSCIJ 


Appears 

Meaning: 


In: POTENT 

MBIAS MBIAS 

E E 

1=1 j=i 


CIJ (I,J) . 
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SUM (3) 


Appears In: COMMON/SUNVEC/ ; INOPT, CHARGE and 

related routines; HIDCEL (see PDIR) . 

Meaning: User-specified x, y and z components, of 

the sun direction vector (direction-of- 
view vector for plots) pointing from 
the satellite center towards the sun 
(or viewer). See DIR(3). SUN values 
will be used to calculate new shadowed 
cell areas only if IAREA > 0. 

SWPCIJ 

Appears In: POTENT 

MBIAS 

Meaning: (-SWPCND(I) + CIJSUM(I)). 

1=1 

SWPSSC 

Appears In: POTENT, OBJSPC 

Meaning: SWPCIJ - SSCIJ. 

T AREA1 

Appears In: HIDCEL 

leaning: Temporary save of AREA1 (IC) for 

printing purposes. 

TIME 

Appears In: CHARGE, TRILIN 

Meaning: The time (in seconds) from the be- 

ginning of cycle 1. 
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UCONP (7) 


Appears In: POTENT and related routines. 

Meaning: Conjugate gradient algorithm u array 

values for each of up to seven con- 
ductors. See subroutine POTENT. 

UDOTAU (double precision) 

Appears In: POTENT and related routines. 

Meaning: Conjugate gradient algorithm dot 

product u n *(Au) n . See subroutine 
POTENT . 

UDOTS 

Appears In: POTENT 

Meaning: Temporary save of first iteration 

UDOTAU for possible convergence 
testing (not implemented) . 

UDO'fU 

Appears In: POTENT, UDOTUC 

Meaning: Conjugate gradient algorithm dot 

product u n *u n . Used in calculation 
of r •? = u n, u n /C* for algorithm one. 

VBIAS (7) 

Appears In: COMMON/MCOND/ ; INKEYW, POTENT, 

PTROUS 

Meaning: If conductor I is biased relative to 

conductor one, then VBIAS (I) is the 
amount of bias. PCOND(I) = PCOND(l) + 
VBIAS (I). Otherwise VBIAS (I) = 0. 
VBIAS is defined in the keyword input 
file. See subroutine INKEYW. 



VFIX (7 ) 


Appears In: COMMON/MCOND/ ; INKEYW, POTENT and 

related routines. 

Meaning: If conductor I is a fixed potential 

conductor then VFIX (I) equals the 
potential at which it is fixed. 
Otherwise VFIX (I) = 0. VFIX is de- 
fined in the keyword input file. 

See subroutine INKEYW. 

vx 

Appears In: HIDCEL, PRSPPJ 

Meaning: QNX2 + 1000. *DIR(1) . x coordinate of 

"end" of direction-of-view vector 
pointing from satellite center to 
viewer. Sec DIR. Used by PRSPPJ 
(entry point PPJSET) . 

VY 

Appears In: HIDCEL, PRSPPJ 

Meaning: QNY2 + DFAC*L>IR(2) . y coordinate of 

"end" of direction-of-view vector. 

See VX. 

VZ 

Appears In: HIDCEL, PRSPPJ 

Meaning: QNZ2 + DFAC"DIR(3) . z coordinate of 

"end" of direction-of-view vector. 

See VZ. 

WPCOND (1024 ) 

- Appears In: OBJDEF, POTENT, TR/LIN 

Meaning: An array containing the potential 

coefficients between surface points 
and their underlying conductors. 
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W4 

Appears In: COPROD, IRCALL, ORCALL and related rou- 

tines. 

Meaning: Geometrical weighting factor determined by 

nesting level of current grid. For grid IG, 
W4 = 2 ** (IG-2) . 

XDIF 

Appear? In: HIDCEL 

Meaning: Span of 2-D x values for current satellite 

plot. Used in defining plot user area. 

XMAX 

Appears In: HIDCEL 

Meaning: Maximum 2-D x value for current satellite 

plot. Used in defining XDIF. 

XMAX1 (600) 

Appears In: HIDCEL and related routines. 

Meaning: XMAXl(I) is the maximum 2-D x value of 

any of the vertices of Al polygon I. 

XMAX2 (150) 

Appears In: HIDCEL and related routines. 

Meaning: XMAX2(I) is the maximum 2-D x value of any 

of the vertices of A2 polygon I. 

XMESH 

Appears In: CHARGE, NASCAP, OBJDEF 

Meaning: The physical size (in meters) of the inner 

grid unit. 
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XMIN 


Appears In: HIDCEL 

Meaning: Minimum 2-D x value for current 

satellite plot. Used in defining 
XDIF. 

XMIN1 (600) 


Appears In: HIDCEL and related routines. 

Meaning: XMINl(I) is the minimum x value of 

any of the vertices of A1 polygon I. 

XMIN2 (150) . _ 

Appears In: HIDCEL and related routines. 

Meaning: XMIN2(I) is the minimum x value of 

any of the vertices of A2 polygon I. 

YCOND ( 7 ) 

Appears In: POTENT and related routines. 

Meaning: Conjugate gradient algorithm y array 

values for each of up to seven con- 
ductors. See subroutine POTENT. 


YDIF 

Appears In: HIDCEL 

Meaning: Span of 2-D y values for current 

satellite plot. Used in defining 
plot user area. 


YMAX 

Appears In: HIDCEL 

Meaning: Maximum 2-D y value for current 

satellite plot. Used in defining 
YDIF. 
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YMAXl (600) 

Appears In: HIDCEL and related routines. 

Meaning: YMAXl (I) is the maximum 2-D y value 

of any of the vertices of A1 polygon I. 

YMAX2 (150) 

Appears In: J.’IDCEL and related routines. 

Meaning: YMAX2 (I) is the maximum 2-D y value 

of any of the vertices of A2 polygon I. 

YMIN 

Appears In: HIDCEL 

Meaning: Minimum 2-D y value for current 

satellite plot. Used in defining 
YD IF. 

YMIN1 (600) 

Appears In: HIDCEL and related routines. 

Meaning: YMINl(I) is the minimum 2-D y value 

of any of the vertices of A1 polygon I. 

YMIN 2 (150) 

Appears In: HIDCEL and related routines. 

Meaning: YMIN2 (I) is the minimum 2-D y value 

of any of the vertices of A2 polygon I. 
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12. SUBROUTINE INDEX 


Subroutine 

Page 

AD DAI 

25Q 

ADDVXB 

233 

APRT 

98 

AREA 

251 

AlCOMP 

252 

AlGEN 

253 

A2GEN 

254 

A2PLOT 

256 

BFIELD 

233 

BKSCAT 

231 

BSCAT 

231 

CALFLX 

222 

CFLICT 

162 

CHARGE 

222 

CMPRSS 

159 

CNVPLT 

106 

CONDUC 

159 

CONPLT 

257 

CONPOT 

258 

COPROD 

189 

CUBE 3 4 

162 

DELETE 

161 

DIAGNO 

159 

DOPLOT 

107 

EFIELD 

234 

ELFLUX 

231 

ELS EC 

231 

ELSEC1 

165 

ESPEC' 

235 

ETNGUN 

225 
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Subroutine 


Page 

FILINP 


235 

FIL111 


163 

FINVER 


259 

FLXDEF 


100 

FREAD 


227 

FSPACE 


227 

GETCEL 


236 

GETNC 


100 

HXDCEL 


261 

IBIT 


167 

ICORN 


191 

INDATA 


100 

INEDGX 


192 

INEDGY 


193 

INEDG2 


194 

INEIMP 


' 165 

INKEYW 

- 

ion 

INOPT 


101 

INPLOT 


107 

INPUT 


161 

INSIDD 


265 

INSIDE 


266 

INSID1 


236 

INTCHK 


267 

INTSEC 


269 

IOFLX 


232 

IPLANX 


195 

IPLANY 


196 

IPLANZ 


197 

IRCALL 


198 

I SPACE 


199 

LCODE 


236 

LINPLN 


270 
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Subroutine 


Page 

LINPLT 


107 

LORDEQ 


161 

LORDER 


101 

MAGCMP 


226 

MATCAL 


271 

MATPLT 


272 

MATPRO 


166 

MOVDAT 


102 

MOVEA1 


273 . 

NASCAP 


96 

NIOOBJ 


163 

NIOOCT 


163 

NIOTET 


163 

NIOWGE 


164 

NUMLTB 


161 

OBJDEF 


162 

OBJSPC 


200 

OCORN 


201 

OCTGON 


164 

OLDEIN 


108 

OPLANX 


202 

OPLANY 


203 

OPLANZ 


204 

ORCALL 


205 

OUTEDG 


206 

PARPLT 


274 

PERMU 


167 

P HOC UR 


238 

PMODQ 


104 

POINT 


108 

POTENT 


207 

PRFLUX 


232 

PROLIN 


275 
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Subroutine 


PROSEC 

PRSPPJ 

PTCOMP 

PTLINE 

PTROUS 

PUPDAT 

PUSH 

QCONCjjL 

QEQN 

QEQN2 

QSPHER 

QSUMER 

RECTAN 

REPLOT 

RES^Q 

REWIND 

ROUSP 

RSDUMP 

RUPDAT 

SATPLT 

JETALL 

SETAN 

SETFAC 

SETFFL 

SETFOC 

SETFTH 

SETFWG 

SETIP 

SETMAX 

SETMC 

SETOP 

SETPOP 

SETWE 


Page 

232 

276 

278 

279 

209 

210 
228 

- 211 . 
166 
166 
164 
104 
164 
108 
104 
103 

280 

103 
212 
281 

104 
228 
282 

285 

286 

287 

288 

105 
230 
162 

105 

106 
230 


Subroutine 

Page 

SHECAL 

239 

SHIELD 

289 

SMFDEF 

230 

SOLFLX 

233 

SOURCE 

226 

SPACE 

213 

SPCPTS 

167 

SPCWGT 

174 

SPECEL 

167 

SRFPTS . , 

174 

SUMQD 

103 

TETRAH 

164 

TRAAN 

237 

TRANA 

236 

TRAND 

237 

TRANE 

237 

TRILIN 

98 

TRNGLS 

165 

UDOTUC 

214 

UDOTUP 

215 

URSETO 

216 

UUPDAT 

217 

UUPDAU 

218 

VADD 

237 

VADDS 

238 

VMULT 

238 

VSET 

238 

WEDGE 

165 

YINIT 

219 

YUPDAT 

220 

ZSYSTM 

16$ 
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